diff --git a/README.md b/README.md
index 3396e8cd..77b68a57 100644
--- a/README.md
+++ b/README.md
@@ -95,7 +95,8 @@ There is a [new documentation website](https://openmopac.github.io) under develo
## Interfaces
While MOPAC is primarily a self-contained command-line program whose behavior is specified by an input file, it also has other modes of
-operation, some of which only require the MOPAC shared library and not the executable.
+operation, some of which only require the MOPAC shared library and not the executable. Note that API calls to the MOPAC library are not
+thread safe. Each thread must load its own instance of the MOPAC library, such as by running independent calling programs.
### MDI Engine
diff --git a/src/MOZYME/iter_for_MOZYME.F90 b/src/MOZYME/iter_for_MOZYME.F90
index 68c7b26e..fffad5e3 100644
--- a/src/MOZYME/iter_for_MOZYME.F90
+++ b/src/MOZYME/iter_for_MOZYME.F90
@@ -15,7 +15,7 @@
! along with this program. If not, see .
subroutine iter_for_MOZYME (ee)
- use molkst_C, only: norbs, step_num, numcal, nscf, escf, &
+ use molkst_C, only: norbs, step_num, numcal, numcal0, nscf, escf, &
& numat, enuclr, atheat, emin, keywrd, moperr, line, use_disk
!
use chanel_C, only: iw, iend, end_fn
@@ -519,7 +519,7 @@ subroutine iter_for_MOZYME (ee)
goto 80
end if
end if
- if (bigscf .or. numcal /= 1) then
+ if (bigscf .or. numcal /= 1+numcal0) then
call diagg (f, nocc1, nvir1, idiagg, partp, indi)
idiagg = idiagg + 1
else
@@ -610,7 +610,7 @@ subroutine iter_for_MOZYME (ee)
backspace (iw)
end if
call isitsc (escf, selcon, emin, iemin, iemax, okscf, niter, itrmax)
- if ( .not. bigscf .and. numcal == 1) then
+ if ( .not. bigscf .and. numcal == 1+numcal0) then
exit
else if (okscf .and. niter > 1 .and. (emin /= 0.d0 .or. niter > 3)) then
exit
diff --git a/src/MOZYME/tidy.F90 b/src/MOZYME/tidy.F90
index 80b215ec..4e892e2a 100644
--- a/src/MOZYME/tidy.F90
+++ b/src/MOZYME/tidy.F90
@@ -16,7 +16,7 @@
subroutine tidy (nmos_loc, nc, ic, n01, c, n02, nnc_loc, ncmo, ln, mn, mode)
use MOZYME_C, only: iorbs, jopt, thresh, numred
- use molkst_C, only: numat, step_num, norbs, moperr, keywrd, numcal, use_disk
+ use molkst_C, only: numat, step_num, step_num0, norbs, moperr, keywrd, numcal, use_disk
use chanel_C, only: iw
implicit none
integer, intent (in) :: mode, n02, nmos_loc
@@ -115,7 +115,7 @@ subroutine tidy (nmos_loc, nc, ic, n01, c, n02, nnc_loc, ncmo, ln, mn, mode)
!
! Do the LMOs of the SCF need to be put at the start of the storage?
!
- if (isnew /= 0 .or. step_num <= 1 .or. step_num == imode(mode)) exit
+ if (isnew /= 0 .or. step_num <= 1+step_num0 .or. step_num == imode(mode)) exit
isnew = 2
if (numred >= numat-1) then
isnew = 1
diff --git a/src/input/getdat.F90 b/src/input/getdat.F90
index 65667df0..8dcea290 100644
--- a/src/input/getdat.F90
+++ b/src/input/getdat.F90
@@ -65,16 +65,6 @@ subroutine getdat(input, output)
call getarg (run, jobnam)
#endif
natoms = 1
- do i = len_trim(jobnam), 1, -1 ! Remove any unprintable characters from the end of the file-name
- if (ichar(jobnam(i:i)) > 39 .and. ichar(jobnam(i:i)) < 126 .or. jobnam(i:i) =="'") exit
- end do
- jobnam(i + 1:) = " "
-!
-! Replace backslash with forward-slash
-!
- do i = 1, len_trim(jobnam)
- if (jobnam(i:i) == backslash) jobnam(i:i) = "/"
- end do
else if (i == 0) then
if (is_PARAM) then
write(line,'(2a)')" PARAM is the parameter optimization program for use with MOPAC"
@@ -99,6 +89,17 @@ subroutine getdat(input, output)
natoms = 1
end if
end if
+! Remove any unprintable characters from the end of the file-name
+ do i = len_trim(jobnam), 1, -1
+ if (ichar(jobnam(i:i)) > 39 .and. ichar(jobnam(i:i)) < 126 .or. jobnam(i:i) =="'") exit
+ end do
+ jobnam(i + 1:) = " "
+!
+! Replace backslash with forward-slash
+!
+ do i = 1, len_trim(jobnam)
+ if (jobnam(i:i) == backslash) jobnam(i:i) = "/"
+ end do
if (natoms == 0) return
!
! Check for the data set in the order: .mop, .dat,
@@ -294,7 +295,7 @@ subroutine getdat(input, output)
i = index(keywrd, "GEO-DAT")
if (i /= 0) keywrd(i:i+6) = "GEO_DAT"
i = index(keywrd, "GEO-REF")
- if (i /= 0) keywrd(i:i+6) = "GEO_REF"
+ if (i /= 0) keywrd(i:i+6) = "GEO_REF"
double_plus = (index(keywrd, " ++ ") /= 0)
if (index(keywrd, " GEO_DAT") + index(keywrd, " SETUP")/= 0) then
nlines = nlines + 3
diff --git a/src/input/getgeo.F90 b/src/input/getgeo.F90
index 1f46c1e6..4d04c0d1 100644
--- a/src/input/getgeo.F90
+++ b/src/input/getgeo.F90
@@ -22,7 +22,7 @@ subroutine getgeo(iread, labels, geo, xyz, lopt, na, nb, nc, int)
use parameters_C, only : ams
!
use molkst_C, only : natoms, keywrd, numat, maxtxt, line, moperr, &
- numcal, id, units, Angstroms, arc_hof_1, arc_hof_2, keywrd_txt, pdb_label
+ numcal, numcal0, id, units, Angstroms, arc_hof_1, arc_hof_2, keywrd_txt, pdb_label
!
use chanel_C, only : iw, ir, input_fn, end_fn, iend
!
@@ -162,11 +162,11 @@ subroutine getgeo(iread, labels, geo, xyz, lopt, na, nb, nc, int)
if (line == '$end') go to 20
if (line(1:1) == '*') go to 20
if (line == ' ') then
- if(natoms == 0 .and. numcal == 1) then
+ if(natoms == 0 .and. numcal == 1+numcal0) then
!
! Check: Is this an ARC file?
!
- numcal = 2
+ numcal = 2+numcal0
rewind (iread)
sum = 0.d0
do i = 1, 10000
@@ -572,7 +572,7 @@ subroutine getgeo(iread, labels, geo, xyz, lopt, na, nb, nc, int)
write(iw,'(/10x,a,i5)')"Faulty atom:", natoms
write(iw,'(/10x,a)')"Faulty line: """//trim(line)//""""
call mopend("Unless MINI is used, optimization flags must be 1, 0, or -1")
- numcal = 2
+ numcal = 2+numcal0
if ((lopt(1,natoms) > 10 .or. lopt(2,natoms) > 10 .or. lopt(3,natoms) > 10) .and. natoms > 1) &
write(iw,'(/10x,a)')" If the geometry is in Gaussian format, add keyword ""AIGIN"" and re-run"
return
@@ -626,7 +626,7 @@ subroutine getgeo(iread, labels, geo, xyz, lopt, na, nb, nc, int)
!***********************************************************************
120 continue
if (natoms == 0) then
- if (numcal == 1) call mopend (' Error detected while reading geometry')
+ if (numcal == 1+numcal0) call mopend (' Error detected while reading geometry')
return
end if
if ( .not. Angstroms) then
diff --git a/src/input/gettxt.F90 b/src/input/gettxt.F90
index db2c2c42..a542e177 100644
--- a/src/input/gettxt.F90
+++ b/src/input/gettxt.F90
@@ -16,7 +16,7 @@
subroutine gettxt
use chanel_C, only: ir, iw, isetup, input_fn
- use molkst_C, only: keywrd, keywrd_quoted, koment, title, refkey, gui, numcal, line, &
+ use molkst_C, only: keywrd, keywrd_quoted, koment, title, refkey, gui, numcal, numcal0, line, &
moperr, allkey, backslash
implicit none
!-----------------------------------------------
@@ -107,7 +107,7 @@ subroutine gettxt
if (.not. exists) then
if (setup_present .and. .not. zero_scf) then
write (line, '(A)') "SETUP FILE """//trim(filen)//""" MISSING."
- numcal = 2
+ numcal = 2+numcal0
if (.not. gui )write(0,'(//30x,a)')' SETUP FILE "'//trim(filen)//'" MISSING'
call mopend (trim(line))
return
@@ -401,7 +401,7 @@ subroutine gettxt
go to 60
50 continue
if (zero_scf) go to 60
- numcal = 2
+ numcal = 2+numcal0
call mopend ('SETUP FILE "'//trim(filen)//'" MISSING')
write(iw,'(a)') " (Setup file name: '"//trim(filen)//"')"
return
@@ -411,7 +411,7 @@ subroutine gettxt
100 continue
call split_keywords(oldkey)
- if (numcal > 1) then
+ if (numcal > 1+numcal0) then
if (index(keywrd,"OLDGEO") /= 0) return ! User forgot to add extra lines for title and comment
if (aux) keywrd = " AUX"
line = "JOB ENDED NORMALLY"
diff --git a/src/input/readmo.F90 b/src/input/readmo.F90
index d51497fd..d5e30203 100644
--- a/src/input/readmo.F90
+++ b/src/input/readmo.F90
@@ -27,7 +27,7 @@ subroutine readmo
!
USE symmetry_C, ONLY: idepfn, locdep, depmul, locpar
!
- use molkst_C, only : ndep, numat, numcal, natoms, nvar, keywrd, dh, &
+ use molkst_C, only : ndep, numat, numcal, numcal0, natoms, nvar, keywrd, dh, &
& verson, is_PARAM, line, nl_atoms, l_feather, backslash, &
& moperr, maxatoms, koment, title, method_pm6, refkey, l_feather_1, &
isok, method_pm6_dh2, caltyp, keywrd_quoted, &
@@ -464,13 +464,13 @@ subroutine readmo
intern = .false.
else
call getgeo (ir, labels, geo, coord, lopt, na, nb, nc, intern)
- if (numcal == 1 .and. natoms == 0) then
+ if (numcal == 1+numcal0 .and. natoms == 0) then
i = index(keywrd, "GEO_DAT")
if (i /= 0) then
write(line,'(2a)')" GEO_DAT file """//trim(line_1)//""" exists, but does not contain any atoms."
write(0,'(//10x,a,//)')trim(line)
call mopend(trim(line))
- else if (.not. gui .and. numcal < 2) then
+ else if (.not. gui .and. numcal < 2+numcal0) then
write(line,'(2a)')" Data set '"//trim(job_fn)//" exists, but does not contain any atoms."
write(0,'(//10x,a,//)')trim(line)
call mopend(trim(line))
@@ -498,7 +498,7 @@ subroutine readmo
coorda(:,:numat) = geo(:,:numat)
numat_old = numat
else if (natoms /= -3) then
- if (moperr .and. numcal == 1) return
+ if (moperr .and. numcal == 1+numcal0) return
if (maxtxt > txtmax) txtmax = maxtxt
txtatm1(:natoms) = txtatm(:natoms)
if (index(keywrd, " RESID") /= 0) txtatm1(:numat)(22:22) = " "
@@ -623,7 +623,7 @@ subroutine readmo
end do
end if
if (natoms < 0 ) then
- if (numcal == 1) rewind ir
+ if (numcal == 1+numcal0) rewind ir
if (.not.isok) then
write (iw, '(A)') &
' Use AIGIN to allow more geometries to be used'
@@ -634,7 +634,7 @@ subroutine readmo
stop
end if
isok = .FALSE.
- if (numcal > 2) then
+ if (numcal > 2+numcal0) then
naigin = naigin + 1
write (iw, '(2/,2A)') ' GAUSSIAN INPUT REQUIRES', &
' STAND-ALONE JOB'
@@ -647,7 +647,7 @@ subroutine readmo
go to 10
end if
end if
- if (natoms == 0 .and. numcal == 1) then
+ if (natoms == 0 .and. numcal == 1+numcal0) then
call mopend ('NO ATOMS IN SYSTEM')
return
end if
@@ -655,7 +655,7 @@ subroutine readmo
!
! Use the old geometry, if one exists
!
- if (numcal == 1) then
+ if (numcal == 1+numcal0) then
write(line,'(a)')" Keyword OLDGEO cannot be used in the first calculation - there is no old geometry"
write(iw,'(//10x,a)')trim(line)
call to_screen(trim(line))
@@ -689,7 +689,7 @@ subroutine readmo
call mopend(trim(line))
return
else
- if (numcal == 1 .and. numat > 50) write(0,'(10x,a)')idate//" Job: '"//trim(jobnam)//"' started successfully"
+ if (numcal == 1+numcal0 .and. numat > 50) write(0,'(10x,a)')idate//" Job: '"//trim(jobnam)//"' started successfully"
end if
end if
maxci = 10000
diff --git a/src/molkst_C.F90 b/src/molkst_C.F90
index 3e5caf45..c5701efa 100644
--- a/src/molkst_C.F90
+++ b/src/molkst_C.F90
@@ -110,6 +110,7 @@ module molkst_C
! Each stage is limited to the same electronic structure. Most calculations
! will only have one stage, e.g. geometry optimization or force constants.
!
+ & job_no0, numcal0, step_num0, & ! Needed for repeated API calls to run_mopac
& mpack, & ! Number of elements in a lower-half-triangle = (norbs*(norbs+1))/2
& n2elec, & ! Number of two-electron integrals
& nscf, & ! Number of SCF calculations done
diff --git a/src/mopend.F90 b/src/mopend.F90
index f21c0515..92eaf314 100644
--- a/src/mopend.F90
+++ b/src/mopend.F90
@@ -42,7 +42,7 @@ subroutine summary(txt, ntxt)
!
!
use chanel_C, only : iw, ir
- use molkst_C, only : line, job_no, natoms, dummy, errtxt
+ use molkst_C, only : line, job_no, job_no0, natoms, dummy, errtxt
implicit none
integer, intent (in) :: ntxt
character, intent (in) :: txt*(*)
@@ -83,7 +83,7 @@ subroutine summary(txt, ntxt)
end if
end if
if (ntxt == 1) then
- if (natoms == 0 .and. job_no == 1) then
+ if (natoms == 0 .and. job_no == job_no0) then
write(iw,'(/10x, a)')"Job failed to run because no atoms were detected in the system"
write(iw,'(10x, a, /)')"The start of the data-set is as follows:"
rewind (ir)
diff --git a/src/output/to_screen.F90 b/src/output/to_screen.F90
index aee4bd63..b9748d30 100644
--- a/src/output/to_screen.F90
+++ b/src/output/to_screen.F90
@@ -63,7 +63,7 @@ subroutine current_version (text)
method_am1, method_mndo, method_pm3, method_rm1, method_mndod, method_pm6, &
method_pm7, nvar, koment, keywrd, zpe, id, density, natoms, formula, press, voigt, &
uhf, nalpha, nbeta, gnorm, mozyme, mol_weight, ilim, &
- line, nscf, time0, sz, ss2, no_pKa, title, jobnam, job_no, fract
+ line, nscf, time0, sz, ss2, no_pKa, title, jobnam, job_no, job_no0, fract
!
use MOZYME_C, only : ncf, ncocc, noccupied, icocc_dim, cocc_dim, nvirtual, icvir_dim, &
nncf, iorbs, cocc, icocc, ncvir, nnce, nce, icvir, cvir, tyres, size_mres, &
@@ -1340,7 +1340,7 @@ subroutine current_version (text)
!
! Don't print processor-independent CPU times for quick jobs - that would waste too much time.
!
- if (time0 > 1.d0 .and. job_no < 4) then
+ if (time0 > 1.d0 .and. job_no < 4+job_no0) then
!
! Deliberately run a time-consuming calculation to work out CPU speed
! The "j" index is set to use up 1.0 seconds on the development computer.
diff --git a/src/output/writmo.F90 b/src/output/writmo.F90
index 9ef3af2b..3421ed70 100644
--- a/src/output/writmo.F90
+++ b/src/output/writmo.F90
@@ -18,7 +18,7 @@ subroutine writmo
use cosmo_C, only : iseps, area, fepsi, cosvol, ediel, solv_energy
!
use molkst_C, only : numat, nclose, nopen, fract, nalpha, nelecs, nbeta, &
- & norbs, nvar, gnorm, iflepo, enuclr,elect, ndep, nscf, numcal, escf, &
+ & norbs, nvar, gnorm, iflepo, enuclr,elect, ndep, nscf, numcal, numcal0, escf, &
& keywrd, os, verson, time0, moperr, last, iscf, id, pressure, mol_weight, &
jobnam, line, mers, uhf, method_indo, &
density, formula, mozyme, mpack, stress, &
@@ -971,7 +971,7 @@ subroutine writmo
call pdbout(31)
close (31)
end if
- if (numcal == 2) then
+ if (numcal == 2+numcal0) then
if (index(keywrd, "OLDGEO") /= 0) then
!
! Write a warning that OLDGEO has been used, so user is aware that multiple ARC files are present
diff --git a/src/run_mopac.F90 b/src/run_mopac.F90
index a3ee681e..0add104e 100644
--- a/src/run_mopac.F90
+++ b/src/run_mopac.F90
@@ -34,7 +34,7 @@ subroutine run_mopac
sparkle, itemp_1, maxtxt, koment, sz, ss2, keywrd_quoted, &
nl_atoms, use_ref_geo, prt_coords, pdb_label, step, &
density, norbs, method_indo, nclose, nopen, backslash, gui, os, git_hash, verson, &
- use_disk, run
+ use_disk, run, numcal0, job_no0, step_num0
!
USE parameters_C, only : tore, ios, iop, iod, eisol, eheat, zs, eheat_sparkles, gss
!
@@ -124,6 +124,10 @@ subroutine run_mopac
stop
endif
end do
+! save state reference to use only relative state values in API calls
+ numcal0 = numcal
+ job_no0 = job_no
+ step_num0 = step_num
!------------------------------------------------------------------------
tore = ios + iop + iod
call fbx ! Factorials and Pascal's triangle (pure constants)
@@ -179,7 +183,7 @@ subroutine run_mopac
10 continue
numcal = numcal + 1 ! A new calculation
job_no = job_no + 1 ! A new job
- if (job_no > 1) then
+ if (job_no > 1+job_no0) then
backspace(ir)
read(ir,'(a)', iostat = i) line
if (i == 0) then
@@ -229,17 +233,17 @@ subroutine run_mopac
l_normal_html = .true.
use_disk = .true.
state_Irred_Rep = " "
- if (job_no > 1) then
+ if (job_no > 1+job_no0) then
i = index(keywrd, " BIGCYCL")
if (i /= 0 .and. index(keywrd,' DRC') == 0) then
i = nint(reada(keywrd, i)) + 1
- if (job_no < i) then
+ if (job_no < i+job_no0) then
fepsi = store_fepsi
goto 90
end if
end if
end if
- if (numcal > 1 .and. numcal < 4 .and. index(keywrd_txt," GEO_DAT") /= 0) then
+ if (numcal > 1+numcal0 .and. numcal < 4+numcal0 .and. index(keywrd_txt," GEO_DAT") /= 0) then
!
! Quickly jump over first three lines
!
@@ -249,7 +253,7 @@ subroutine run_mopac
natoms = i
call gettxt
end if
- if (numcal > 1) call to_screen("To_file: Leaving MOPAC")
+ if (numcal > 1+numcal0) call to_screen("To_file: Leaving MOPAC")
!
! Read in all the data for the current job
!
@@ -265,12 +269,12 @@ subroutine run_mopac
if (j /= 0) i = i - 6 + j
end if
inquire(file=line(:i)//".den", exist=l_OLDDEN)
-90 if (moperr .and. numcal == 1 .and. natoms > 1) goto 101
- if (moperr .and. numcal == 1 .and. index(keywrd_txt," GEO_DAT") == 0) goto 100
+90 if (moperr .and. numcal == 1+numcal0 .and. natoms > 1) goto 101
+ if (moperr .and. numcal == 1+numcal0 .and. index(keywrd_txt," GEO_DAT") == 0) goto 100
if (moperr) goto 101
! Adjust maximum number of threads using the OpenMP API
#ifdef _OPENMP
- if (numcal == 1) default_num_threads = omp_get_max_threads()
+ if (numcal == 1+numcal0) default_num_threads = omp_get_max_threads()
i = index(keywrd, " THREADS")
if (i > 0) then
num_threads = nint(reada(keywrd, i))
@@ -279,7 +283,7 @@ subroutine run_mopac
end if
call omp_set_num_threads(num_threads)
#endif
- if (numcal == 1) then
+ if (numcal == 1+numcal0) then
#ifdef MKL
num_threads = min(mkl_get_max_threads(), 20)
i = index(keywrd, " THREADS")
@@ -357,7 +361,7 @@ subroutine run_mopac
lgpu = (lgpu_ref .and. natoms > 100) ! Warning - there are problems with UHF calculations on small systems
#endif
end if
- if (.not. gui .and. numcal == 1 .and. natoms == 0) then
+ if (.not. gui .and. numcal == 1+numcal0 .and. natoms == 0) then
write(line,'(2a)')" Data set exists, but does not contain any atoms."
write(0,'(//10x,a,//)')trim(line)
call mopend(trim(line))
@@ -376,7 +380,7 @@ subroutine run_mopac
natoms = 0
goto 100
end if
- if (numcal == 1 .and. moperr .or. natoms == 0) then
+ if (numcal == 1+numcal0 .and. moperr .or. natoms == 0) then
!
! Check for spurious "extra" data
!