Skip to content

Commit

Permalink
Initial commit of OCv1a from OpenCalphad website.
Browse files Browse the repository at this point in the history
  • Loading branch information
richardotis committed Aug 20, 2014
1 parent cb48bd9 commit dff789d
Show file tree
Hide file tree
Showing 53 changed files with 48,488 additions and 0 deletions.
Empty file added .gitignore
Empty file.
32 changes: 32 additions & 0 deletions Makefile
Original file line number Diff line number Diff line change
@@ -0,0 +1,32 @@
OBJS=metlib3.o tpfun4.o lukasnum.o pmod25.o matsmin.o smp1.o pmon6.o
EXE=oc1A

all:
gfortran -o linkoc linkocdate.F90
./linkoc
make $(EXE)


metlib3.o: utilities/metlib3.F90
gfortran -c -fbounds-check utilities/metlib3.F90

tpfun4.o: utilities/tpfun4.F90
gfortran -c -fbounds-check utilities/tpfun4.F90

lukasnum.o: numlib/lukasnum.F90
gfortran -c -fbounds-check numlib/lukasnum.F90

pmod25.o: models/pmod25.F90
gfortran -c -fbounds-check models/pmod25.F90

matsmin.o: minimizer/matsmin.F90
gfortran -c -fbounds-check minimizer/matsmin.F90

smp1.o: stepmapplot/smp1.F90
gfortran -c -fbounds-check stepmapplot/smp1.F90

pmon6.o: userif/pmon6.F90
gfortran -c -fbounds-check userif/pmon6.F90

$(EXE): $(OBJS) linkocdate.F90
gfortran -o $(EXE) -fbounds-check pmain1.F90 $(OBJS)
178 changes: 178 additions & 0 deletions TQlib/crfe.TDB
Original file line number Diff line number Diff line change
@@ -0,0 +1,178 @@

$ Database file written 2012- 9- 7
$ From database: SSOL2
ELEMENT /- ELECTRON_GAS 0.0000E+00 0.0000E+00 0.0000E+00!
ELEMENT VA VACUUM 0.0000E+00 0.0000E+00 0.0000E+00!
ELEMENT CR BCC_A2 5.1996E+01 4.0500E+03 2.3560E+01!
ELEMENT FE BCC_A2 5.5847E+01 4.4890E+03 2.7280E+01!


FUNCTION GHSERCR 2.98150E+02 -8856.94+157.48*T-26.908*T*LN(T)
+.00189435*T**2-1.47721E-06*T**3+139250*T**(-1); 2.18000E+03 Y
-34869.344+344.18*T-50*T*LN(T)-2.88526E+32*T**(-9); 6.00000E+03 N !
FUNCTION GPCRLIQ 2.98150E+02 +YCRLIQ#*EXP(ZCRLIQ#); 6.00000E+03 N !
FUNCTION GFELIQ 2.98150E+02 +12040.17-6.55843*T-3.6751551E-21*T**7
+GHSERFE#; 1.81100E+03 Y
-10839.7+291.302*T-46*T*LN(T); 6.00000E+03 N !
FUNCTION GPFELIQ 2.98150E+02 +YFELIQ#*EXP(ZFELIQ#); 6.00000E+03 N !
FUNCTION GPCRBCC 2.98150E+02 +YCRBCC#*EXP(ZCRBCC#); 6.00000E+03 N !
FUNCTION GHSERFE 2.98150E+02 +1225.7+124.134*T-23.5143*T*LN(T)
-.00439752*T**2-5.8927E-08*T**3+77359*T**(-1); 1.81100E+03 Y
-25383.581+299.31255*T-46*T*LN(T)+2.29603E+31*T**(-9); 6.00000E+03 N
!
FUNCTION GPFEBCC 2.98150E+02 +YFEBCC#*EXP(ZFEBCC#); 6.00000E+03 N !
FUNCTION GCRFCC 2.98150E+02 +7284+.163*T+GHSERCR#; 6.00000E+03 N !
FUNCTION GFEFCC 2.98150E+02 -1462.4+8.282*T-1.15*T*LN(T)+6.4E-04*T**2
+GHSERFE#; 1.81100E+03 Y
-27098.266+300.25256*T-46*T*LN(T)+2.78854E+31*T**(-9); 6.00000E+03 N
!
FUNCTION GPFEFCC 2.98150E+02 +YFEFCC#*EXP(ZFEFCC#); 6.00000E+03 N !
FUNCTION GPSIG1 2.98150E+02 +1.09E-04*P; 6.00000E+03 N !
FUNCTION GPSIG2 2.98150E+02 +1.117E-04*P; 6.00000E+03 N !
FUNCTION YCRLIQ 2.98150E+02 +VCRLIQ#*EXP(-ECRLIQ#); 6.00000E+03 N !
FUNCTION ZCRLIQ 2.98150E+02 +1*LN(XCRLIQ#); 6.00000E+03 N !
FUNCTION YFELIQ 2.98150E+02 +VFELIQ#*EXP(-EFELIQ#); 6.00000E+03 N !
FUNCTION ZFELIQ 2.98150E+02 +1*LN(XFELIQ#); 6.00000E+03 N !
FUNCTION YCRBCC 2.98150E+02 +VCRBCC#*EXP(-ECRBCC#); 6.00000E+03 N !
FUNCTION ZCRBCC 2.98150E+02 +1*LN(XCRBCC#); 6.00000E+03 N !
FUNCTION YFEBCC 2.98150E+02 +VFEBCC#*EXP(-EFEBCC#); 6.00000E+03 N !
FUNCTION ZFEBCC 2.98150E+02 +1*LN(XFEBCC#); 6.00000E+03 N !
FUNCTION YFEFCC 2.98150E+02 +VFEFCC#*EXP(-EFEFCC#); 6.00000E+03 N !
FUNCTION ZFEFCC 2.98150E+02 +1*LN(XFEFCC#); 6.00000E+03 N !
FUNCTION VCRLIQ 2.98150E+02 +7.653E-06*EXP(ACRLIQ#); 6.00000E+03 N
!
FUNCTION ECRLIQ 2.98150E+02 +1*LN(CCRLIQ#); 6.00000E+03 N !
FUNCTION XCRLIQ 2.98150E+02 +1*EXP(.8*DCRLIQ#)-1; 6.00000E+03 N !
FUNCTION VFELIQ 2.98150E+02 +6.46677E-06*EXP(AFELIQ#); 6.00000E+03
N !
FUNCTION EFELIQ 2.98150E+02 +1*LN(CFELIQ#); 6.00000E+03 N !
FUNCTION XFELIQ 2.98150E+02 +1*EXP(.8484467*DFELIQ#)-1; 6.00000E+03
N !
FUNCTION VCRBCC 2.98150E+02 +7.188E-06*EXP(ACRBCC#); 6.00000E+03 N
!
FUNCTION ECRBCC 2.98150E+02 +1*LN(CCRBCC#); 6.00000E+03 N !
FUNCTION XCRBCC 2.98150E+02 +1*EXP(.8*DCRBCC#)-1; 6.00000E+03 N !
FUNCTION VFEBCC 2.98150E+02 +7.042095E-06*EXP(AFEBCC#); 6.00000E+03
N !
FUNCTION EFEBCC 2.98150E+02 +1*LN(CFEBCC#); 6.00000E+03 N !
FUNCTION XFEBCC 2.98150E+02 +1*EXP(.7874195*DFEBCC#)-1; 6.00000E+03
N !
FUNCTION VFEFCC 2.98150E+02 +6.688726E-06*EXP(AFEFCC#); 6.00000E+03
N !
FUNCTION EFEFCC 2.98150E+02 +1*LN(CFEFCC#); 6.00000E+03 N !
FUNCTION XFEFCC 2.98150E+02 +1*EXP(.8064454*DFEFCC#)-1; 6.00000E+03
N !
FUNCTION ACRLIQ 2.98150E+02 +1.7E-05*T+9.2E-09*T**2; 6.00000E+03 N
!
FUNCTION CCRLIQ 2.98150E+02 3.72E-11; 6.00000E+03 N !
FUNCTION DCRLIQ 2.98150E+02 +1*LN(BCRLIQ#); 6.00000E+03 N !
FUNCTION AFELIQ 2.98150E+02 +1.135E-04*T; 6.00000E+03 N !
FUNCTION CFELIQ 2.98150E+02 +4.22534787E-12+2.71569924E-14*T;
6.00000E+03 N !
FUNCTION DFELIQ 2.98150E+02 +1*LN(BFELIQ#); 6.00000E+03 N !
FUNCTION ACRBCC 2.98150E+02 +1.7E-05*T+9.2E-09*T**2; 6.00000E+03 N
!
FUNCTION CCRBCC 2.98150E+02 2.08E-11; 6.00000E+03 N !
FUNCTION DCRBCC 2.98150E+02 +1*LN(BCRBCC#); 6.00000E+03 N !
FUNCTION AFEBCC 2.98150E+02 +2.3987E-05*T+1.2845E-08*T**2;
6.00000E+03 N !
FUNCTION CFEBCC 2.98150E+02 +2.20949565E-11+2.41329523E-16*T;
6.00000E+03 N !
FUNCTION DFEBCC 2.98150E+02 +1*LN(BFEBCC#); 6.00000E+03 N !
FUNCTION AFEFCC 2.98150E+02 +7.3097E-05*T; 6.00000E+03 N !
FUNCTION CFEFCC 2.98150E+02 +2.62285341E-11+2.71455808E-16*T;
6.00000E+03 N !
FUNCTION DFEFCC 2.98150E+02 +1*LN(BFEFCC#); 6.00000E+03 N !
FUNCTION BCRLIQ 2.98150E+02 +1+4.65E-11*P; 6.00000E+03 N !
FUNCTION BFELIQ 2.98150E+02 +1+4.98009787E-12*P+3.20078924E-14*T*P;
6.00000E+03 N !
FUNCTION BCRBCC 2.98150E+02 +1+2.6E-11*P; 6.00000E+03 N !
FUNCTION BFEBCC 2.98150E+02 +1+2.80599565E-11*P+3.06481523E-16*T*P;
6.00000E+03 N !
FUNCTION BFEFCC 2.98150E+02 +1+3.25236341E-11*P+3.36607808E-16*T*P;
6.00000E+03 N !
FUNCTION UN_ASS 298.15 0; 300 N !

TYPE_DEFINITION % SEQ *!
DEFINE_SYSTEM_DEFAULT ELEMENT 2 !
DEFAULT_COMMAND DEF_SYS_ELEMENT VA /- !


PHASE LIQUID:L % 1 1.0 !
CONSTITUENT LIQUID:L :CR,FE : !

PARAMETER G(LIQUID,CR;0) 2.98150E+02 +24339.955-11.420225*T
+2.37615E-21*T**7+GHSERCR#+GPCRLIQ#; 2.18000E+03 Y
+18409.36-8.563683*T+2.88526E+32*T**(-9)+GHSERCR#+GPCRLIQ#; 6.00000E+03
N REF283 !
PARAMETER G(LIQUID,FE;0) 2.98150E+02 +GFELIQ#+GPFELIQ#; 6.00000E+03
N REF283 !
PARAMETER G(LIQUID,CR,FE;0) 2.98150E+02 -14550+6.65*T; 6.00000E+03
N REF107 !


TYPE_DEFINITION & GES A_P_D BCC_A2 MAGNETIC -1.0 4.00000E-01 !
PHASE BCC_A2 %& 2 1 3 !
CONSTITUENT BCC_A2 :CR%,FE% : VA% : !

PARAMETER G(BCC_A2,CR:VA;0) 2.98150E+02 +GHSERCR#+GPCRBCC#;
6.00000E+03 N REF283 !
PARAMETER TC(BCC_A2,CR:VA;0) 2.98150E+02 -311.5; 6.00000E+03 N
REF281 !
PARAMETER BMAGN(BCC_A2,CR:VA;0) 2.98150E+02 -.01; 6.00000E+03 N
REF281 !
PARAMETER G(BCC_A2,FE:VA;0) 2.98150E+02 +GHSERFE#+GPFEBCC#;
6.00000E+03 N REF283 !
PARAMETER TC(BCC_A2,FE:VA;0) 2.98150E+02 1043; 6.00000E+03 N REF281 !
PARAMETER BMAGN(BCC_A2,FE:VA;0) 2.98150E+02 2.22; 6.00000E+03 N
REF281 !
PARAMETER G(BCC_A2,CR,FE:VA;0) 2.98150E+02 +20500-9.68*T; 6.00000E+03
N REF107 !
PARAMETER TC(BCC_A2,CR,FE:VA;0) 2.98150E+02 1650; 6.00000E+03 N
REF107 !
PARAMETER TC(BCC_A2,CR,FE:VA;1) 2.98150E+02 550; 6.00000E+03 N
REF107 !
PARAMETER BMAGN(BCC_A2,CR,FE:VA;0) 2.98150E+02 -.85; 6.00000E+03 N
REF107 !


TYPE_DEFINITION ' GES A_P_D FCC_A1 MAGNETIC -3.0 2.80000E-01 !
PHASE FCC_A1 %' 2 1 1 !
CONSTITUENT FCC_A1 :CR,FE% : VA% : !

PARAMETER G(FCC_A1,CR:VA;0) 2.98150E+02 +GCRFCC#+GPCRBCC#;
6.00000E+03 N REF281 !
PARAMETER TC(FCC_A1,CR:VA;0) 2.98150E+02 -1109; 6.00000E+03 N
REF281 !
PARAMETER BMAGN(FCC_A1,CR:VA;0) 2.98150E+02 -2.46; 6.00000E+03 N
REF281 !
PARAMETER G(FCC_A1,FE:VA;0) 2.98150E+02 +GFEFCC#+GPFEFCC#;
6.00000E+03 N REF283 !
PARAMETER TC(FCC_A1,FE:VA;0) 2.98150E+02 -201; 6.00000E+03 N REF281 !
PARAMETER BMAGN(FCC_A1,FE:VA;0) 2.98150E+02 -2.1; 6.00000E+03 N
REF281 !
PARAMETER G(FCC_A1,CR,FE:VA;0) 2.98150E+02 +10833-7.477*T;
6.00000E+03 N REF107 !
PARAMETER G(FCC_A1,CR,FE:VA;1) 2.98150E+02 1410; 6.00000E+03 N
REF107 !


PHASE SIGMA % 3 8 4 18 !
CONSTITUENT SIGMA :FE : CR : CR,FE : !

PARAMETER G(SIGMA,FE:CR:CR;0) 2.98150E+02 +8*GFEFCC#+22*GHSERCR#+92300
-95.96*T+GPSIG1#; 6.00000E+03 N REF107 !
PARAMETER G(SIGMA,FE:CR:FE;0) 2.98150E+02 +8*GFEFCC#+4*GHSERCR#
+18*GHSERFE#+117300-95.96*T+GPSIG2#; 6.00000E+03 N REF107 !

LIST_OF_REFERENCES
NUMBER SOURCE
REF283 'Alan Dinsdale, SGTE Data for Pure Elements,
Calphad Vol 15(1991) p 317-425,
also in NPL Report DMA(A)195 Rev. August 1990'
REF281 'Alan Dinsdale, SGTE Data for Pure Elements, NPL Report DMA(A)195
September 1989'
REF107 'J-O Andersson, B. Sundman, CALPHAD Vol 11, (1987), p 83-92
TRITA 0270 (1986); CR-FE'
!

2 changes: 2 additions & 0 deletions TQlib/linktest.txt
Original file line number Diff line number Diff line change
@@ -0,0 +1,2 @@
gfortran -o test1 testtq.F90 tqlib.o oclib.a

158 changes: 158 additions & 0 deletions TQlib/testtq.F90
Original file line number Diff line number Diff line change
@@ -0,0 +1,158 @@
! first test program of OC-TQ
program octq1
use octq
implicit none
! maxel and maxph defined in pmod package
! integer, parameter :: maxel=10,maxph=20
integer n,nel,n1,n2,n3,n4,nph,ics,ip,cnum(maxel+3),mm,m2
character filename*60,elnames(maxel)*2,phnames(maxph)*24
character condition*60,line*80,statevar*60,quest*60,ch1*1
character target*60
double precision value,dummy,tp(2),mel(maxel)
double precision xf(maxel),pxf(10*maxph),npf(maxph)
type(gtp_equilibrium_data), pointer :: ceq
! set some defaults
n=0
filename='crfe '
! initiate
call tqini(n,ceq)
if(gx%bmperr.ne.0) goto 1000
! read database file
call tqrfil(filename,ceq)
if(gx%bmperr.ne.0) goto 1000
! find out about the system:
! number of elements and their names
call tqgcom(nel,elnames,ceq)
if(gx%bmperr.ne.0) goto 1000
! number of phases and their names
call tqgnp(nph,ceq)
if(gx%bmperr.ne.0) goto 1000
do n=1,nph
call tqgpn(n,phnames(n),ceq)
if(gx%bmperr.ne.0) goto 1000
enddo
! -------------------------------------
write(*,10)nel,(elnames(n),n=1,nel)
10 format(/'System with ',i2,' elements: ',10(a,', '))
write(*,20)nph,(phnames(n)(1:len_trim(phnames(n))),n=1,nph)
20 format('and ',i3,' phases: ',10(a,', '))
! --------------------------------------
! set default values
tp(1)=1.0D3
tp(2)=1.0D5
do n=1,nel
xf(n)=1.0/dble(nel)
enddo
! ask for conditions
100 continue
write(*,105)
105 format(/'Give conditions:')
ip=len(line)
dummy=tp(1)
call gparrd('Temperature: ',line,ip,tp(1),dummy,nohelp)
if(buperr.ne.0) goto 1000
if(tp(1).lt.1.0d0) then
write(*,*)'Temperature must be larger than 1 K'
tp(1)=1.0D0
endif
dummy=tp(2)
call gparrd('Pressure: ',line,ip,tp(2),dummy,nohelp)
if(buperr.ne.0) goto 1000
if(tp(2).lt.1.0d0) then
write(*,*)'Pressure must be larger than 1 Pa'
tp(2)=1.0D0
endif
do n=1,nel-1
quest='Mole fraction of '//elnames(n)(1:len_trim(elnames(n)))//':'
dummy=xf(n)
call gparrd(quest,line,ip,xf(n),dummy,nohelp)
if(buperr.ne.0) goto 1000
if(xf(n).lt.1.0d-6) then
write(*,*)'Fraction must be larger than 1.0D-6'
xf(n)=1.0D-6
endif
enddo
! -------------------------------------
! set conditions
n1=0
n2=0
condition='T'
call tqsetc(condition,n1,n2,tp(1),cnum(1),ceq)
if(gx%bmperr.ne.0) goto 1000
condition='P'
call tqsetc(condition,n1,n2,tp(2),cnum(2),ceq)
if(gx%bmperr.ne.0) goto 1000
condition='N'
call tqsetc(condition,n1,n2,one,cnum(3),ceq)
if(gx%bmperr.ne.0) goto 1000
do n=1,nel-1
condition='X'
call tqsetc(condition,n,n2,xf(n),cnum(3+n),ceq)
if(gx%bmperr.ne.0) goto 1000
enddo
!--------------------------------------
! calculate the equilibria
target=' '
n1=0
n2=0
call tqce(target,n1,n2,value,ceq)
if(gx%bmperr.ne.0) goto 1000
write(*,300)
300 format(/'Successful calculation')
!--------------------------------------
! list some results
! amount of all phases
statevar='NP'
n1=-1
n2=0
n3=size(npf)
call tqgetv(statevar,n1,n2,n3,npf,ceq)
if(gx%bmperr.ne.0) goto 1000
write(*,505)n3,(npf(n),n=1,n3)
505 format(/'Amount of ',i2,' phases: ',(10F7.4))
mm=0
phloop: do n=1,nph
icsloop: do ics=1,noofcs(n)
mm=mm+1
if(npf(mm).gt.zero) then
! the phase is stable if it has a positive amount
if(ics.gt.1) then
write(*,510)phnames(n)(1:len_trim(phnames(n))),ics,npf(mm)
else
write(*,511)phnames(n)(1:len_trim(phnames(n))),npf(mm)
endif
510 format(/'Stable phase: ',a,'#',i1,', amount: ',1PE12.4,&
', mole fractions:',3i3)
511 format(/'Stable phase: ',a,', amount: ',1PE12.4,&
', mole fractions:',3i3)
! composition of stable phase, n2=-1 means all fractions
statevar='X'
n2=-1
! Use extended phase index: 10*phase number + compset number ...
n4=size(pxf)
call tqgetv(statevar,10*n+ics,n2,n4,pxf,ceq)
write(*,520)(elnames(m2),pxf(m2),m2=1,n4)
520 format(3(a,': ',F9.6,', '))
endif
enddo icsloop
enddo phloop
!--------------------------------------
! loop
write(*,*)
ip=len(line)
call gparcd('Any more calculations?',line,ip,1,ch1,'N',nohelp)
if(ch1.ne.'N') goto 100
!--------------------------------------
1000 continue
if(gx%bmperr.ne.0) then
if(gx%bmperr.ge.4000 .and. gx%bmperr.le.4220) then
write(*,1010)gx%bmperr,bmperrmess(gx%bmperr)
1010 format(' *** Error ',i5/a)
else
write(*,1020)gx%bmperr
1020 format(' *** Error ',i5/'Unknown reason')
endif
endif
write(*,*)
write(*,*)'Auf wiedersehen'
end program octq1
Loading

0 comments on commit dff789d

Please sign in to comment.