c program to read the original thermodynamic datafile of H & P 1990 program trans open (9, file='pow91', status= 'old') open (10,file='new.dat') call rdin end subroutine rdin c this subroutine reads part of the H&P data character text(132)*1, name*8, cnum*8, gnum*12, twod*2 character snum*8, vnum*8 real rnum, comp(12), rgib, g, reas, s, reav real v, a, b, c, e, b2, b4,l1, l2 integer index, igst, igen, isst, isen, ivst, iven integer ihptov(12) data ihptov/4,7,3,9,2,8,6,1,5,10,11,12/ 1000 format(132a1) 2000 format (a8,4(i2),' H= ',g13.7) 2002 format (5(g13.7,1x)) 2001 format (12(f5.2,1x)) c entropies of the elements at 298, 1bar sk = 64.68 sna = 51.30 sca = 41.63 sc = 5.74 sti = 30.63 sal = 28.35 ssi = 18.81 smg = 32.68 smn = 32.01 sfe = 27.28 c values from Robie so = 205.20 /2. sh = 130.70 /2. 5 read(9,1000,end=900) text c inst and inen specify the positions of initial and terminal * character of the phase name in the array called text * the phase name is then written to the variable name inst = 0 inen = 0 do i = 1, 132 if (text(i) .ne. ' ') then inst = i goto 10 end if end do 10 do i = inst, 132 if (text(i) .eq. ' ') then inen = i-1 goto 20 end if end do 20 write(name,1000) (text(i), i = inst, inen) c reading indices and real volues of components from the array called text * indices are read and then written into the variable index as integers. * Then the value for the component corresponding to the actual index * is read and written ("internal write" statement) to the varaible cnum. * cnum is written to rnum which is then assigned to the corresponding * element of the array compon. do i = 1, 12 nst = 15 + (i-1)*9 if (text(nst).eq.'0'.and.text(nst-1).eq.' ') goto 50 write(twod,1000) (text(j), j = nst-1,nst) read(twod,*) index write(cnum,1000) (text(j), j = nst+1,nst+7) read(cnum,*) rnum comp(ihptov(index)) = rnum end do c finding, reading and writing the values for g 50 do i=nst+1, 132 if(text(i).ne.' ') then igst = i goto 60 end if end do 60 do i=igst+1, 132 if(text(i).eq.' ') then igen = i-1 goto 70 end if end do 70 write(gnum,1000)(text(j),j=igst, igen) read(gnum,*) rgib g = rgib c finding, reading and writing values for s do i = igen+1, 132 if(text(i).ne.' ') then isst = i goto 80 end if end do 80 do i = isst+1, 132 if(text(i).eq.' ') then isen = i-1 goto 90 end if end do 90 write(snum,1000)(text(j), j = isst, isen) read(snum,*) reas s = reas c finding, reading and writing the values for v do i = isen+1, 132 if(text(i).ne.' ') then ivst = i goto 100 end if end do 100 do i = ivst+1, 132 if(text(i).eq.' ') then iven = i-1 goto 110 end if end do 110 write(vnum,1000)(text(j), j = ivst, iven) read(vnum,*) reav v = reav c end of first line in h&p data file c as the second line of the data file has a fixed format, it will be read in a standard c fortran way read(9,*) a, b, c, e, b2, b4, lam if(lam.ne.0) then backspace(9) read(9,*) a, b, c, e, b2, b4, lam, l1, l2 end if ilam = lam * 10 g = g * 1000. s = s * 1000. a = a * 1000. b = b * 1000. c = c * 1000. d = d * 1000. e = e * 1000. b4 = -b4 / 1000. write (10,2000) name,1,0,ilam,0,g c convert hp H to G dsf = s - comp(1) * sna - comp(2) * smg * - comp(3) * sal - comp(4) * ssi * - comp(5) * sk - comp(6) * sca * - comp(7) * sti - comp(8) * smn * - comp(9) * sfe - comp(10) * so * - comp(11) * sh - comp(12) * sc g = g - 298.15 * dsf c convert to oxide stoichiometry comp(1) = comp(1)/2. comp(3) = comp(3)/2. comp(5) = comp(5)/2. comp(11) = comp(11)/2. comp(10) = comp(10) - comp(1) - comp(2) - 3. * comp(3) * - 2. * comp(4) - comp(5) - comp(6) * - 2. * comp(7) - comp(8) - comp(9) * - comp(11) - 2. * comp(12) comp(10) = comp(10) / 2. write (10,2001) (comp(j), j = 1, 12) write (10,2002) g,s,v,a,b,c,0.,e,0.,0.,b2,0.,b4,0.,0.,0. if (ilam.ne.0) then write (10,2002) l1,l2*1000.,0.,0.,0.,0.,0.,0.,0. end if do i = 1, 12 comp(i) = 0. end do goto 5 900 end