%!PS-Adobe-3.0 %%Title: tmpa %%For: David Jeffery %%Creator: a2ps version 4.13 %%CreationDate: Tue Jun 12 11:57:55 2007 %%BoundingBox: 24 24 588 768 %%DocumentData: Clean7Bit %%Orientation: Landscape %%Pages: 7 %%PageOrder: Ascend %%DocumentMedia: Letter 612 792 0 () () %%DocumentNeededResources: font Courier %%+ font Courier-Bold %%+ font Courier-BoldOblique %%+ font Courier-Oblique %%+ font Helvetica %%+ font Helvetica-Bold %%+ font Symbol %%+ font Times-Bold %%+ font Times-Roman %%DocumentProcessColors: Black %%DocumentSuppliedResources: procset a2ps-a2ps-hdr %%+ procset a2ps-black+white-Prolog %%+ encoding ISO-8859-1Encoding %%EndComments /a2psdict 200 dict def a2psdict begin %%BeginProlog %%Copyright: (c) 1988, 89, 90, 91, 92, 93 Miguel Santana %%Copyright: (c) 1995, 96, 97, 98 Akim Demaille, Miguel Santana % Check PostScript language level. /languagelevel where { pop /gs_languagelevel languagelevel def } { /gs_languagelevel 1 def } ifelse % EPSF import as in the Red Book /BeginInclude { /b4_Inc_state save def % Save state for cleanup /dict_count countdictstack def % Count objects on dict stack /op_count count 1 sub def % Count objects on operand stack userdict begin 0 setgray 0 setlinecap 1 setlinewidth 0 setlinejoin 10 setmiterlimit [ ] 0 setdash newpath gs_languagelevel 1 ne { false setstrokeadjust false setoverprint } if } bind def /EndInclude { count op_count sub { pos } repeat % Clean up stacks countdictstack dict_count sub { end } repeat b4_Inc_state restore } bind def /BeginEPSF { BeginInclude /showpage { } def } bind def /EndEPSF { EndInclude } bind def % Page prefeed /page_prefeed { % bool -> - statusdict /prefeed known { statusdict exch /prefeed exch put } { pop } ifelse } bind def /deffont { findfont exch scalefont def } bind def /reencode_font { findfont reencode 2 copy definefont pop def } bind def % Function c-show (str => -) % centers text only according to x axis. /c-show { dup stringwidth pop 2 div neg 0 rmoveto show } bind def % Function l-show (str => -) % prints texts so that it ends at currentpoint /l-show { dup stringwidth pop neg 0 rmoveto show } bind def % center-fit show (str w => -) % show centered, and scale currentfont so that the width is less than w /cfshow { exch dup stringwidth pop % If the title is too big, try to make it smaller 3 2 roll 2 copy gt { % if, i.e. too big exch div currentfont exch scalefont setfont } { % ifelse pop pop } ifelse c-show % center title } bind def % Return the y size of the current font % - => fontsize /currentfontsize { currentfont /FontMatrix get 3 get 1000 mul } bind def % reencode the font % -> /reencode { %def dup length 5 add dict begin { %forall 1 index /FID ne { def }{ pop pop } ifelse } forall /Encoding exch def % Use the font's bounding box to determine the ascent, descent, % and overall height; don't forget that these values have to be % transformed using the font's matrix. % We use `load' because sometimes BBox is executable, sometimes not. % Since we need 4 numbers an not an array avoid BBox from being executed /FontBBox load aload pop FontMatrix transform /Ascent exch def pop FontMatrix transform /Descent exch def pop /FontHeight Ascent Descent sub def % Define these in case they're not in the FontInfo (also, here % they're easier to get to. /UnderlinePosition 1 def /UnderlineThickness 1 def % Get the underline position and thickness if they're defined. currentdict /FontInfo known { FontInfo dup /UnderlinePosition known { dup /UnderlinePosition get 0 exch FontMatrix transform exch pop /UnderlinePosition exch def } if dup /UnderlineThickness known { /UnderlineThickness get 0 exch FontMatrix transform exch pop /UnderlineThickness exch def } if } if currentdict end } bind def % Function print line number ( # -) /# { gsave sx cw mul neg 2 div 0 rmoveto f# setfont c-show grestore } bind def % -------- Some routines to enlight plain b/w printings --------- % Underline % width -- /dounderline { currentpoint gsave moveto 0 currentfont /Descent get currentfontsize mul rmoveto 0 rlineto stroke grestore } bind def % Underline a string % string -- /dounderlinestring { stringwidth pop dounderline } bind def /UL { /ul exch store } bind def % Draw a box of WIDTH wrt current font % width -- /dobox { currentpoint gsave newpath moveto 0 currentfont /Descent get currentfontsize mul rmoveto dup 0 rlineto 0 currentfont /FontHeight get currentfontsize mul rlineto neg 0 rlineto closepath stroke grestore } bind def /BX { /bx exch store } bind def % Box a string % string -- /doboxstring { stringwidth pop dobox } bind def % % ------------- Color routines --------------- % /FG /setrgbcolor load def % Draw the background % width -- /dobackground { currentpoint gsave newpath moveto 0 currentfont /Descent get currentfontsize mul rmoveto dup 0 rlineto 0 currentfont /FontHeight get currentfontsize mul rlineto neg 0 rlineto closepath bgcolor aload pop setrgbcolor fill grestore } bind def % Draw bg for a string % string -- /dobackgroundstring { stringwidth pop dobackground } bind def /BG { dup /bg exch store { mark 4 1 roll ] /bgcolor exch store } if } bind def /Show { bg { dup dobackgroundstring } if ul { dup dounderlinestring } if bx { dup doboxstring } if show } bind def % Function T(ab), jumps to the n-th tabulation in the current line /T { cw mul x0 add bg { dup currentpoint pop sub dobackground } if ul { dup currentpoint pop sub dounderline } if bx { dup currentpoint pop sub dobox } if y0 moveto } bind def % Function n: move to the next line /n { /y0 y0 bfs sub store x0 y0 moveto } bind def % Function N: show and move to the next line /N { Show /y0 y0 bfs sub store x0 y0 moveto } bind def /S { Show } bind def %%BeginResource: procset a2ps-a2ps-hdr 2.0 2 %%Copyright: (c) 1988, 89, 90, 91, 92, 93 Miguel Santana %%Copyright: (c) 1995, 96, 97, 98 Akim Demaille, Miguel Santana % Function title: prints page header. % are passed as argument /title { % 1. Draw the background x v get y v get moveto gsave 0 th 2 div neg rmoveto th setlinewidth 0.95 setgray pw 0 rlineto stroke grestore % 2. Border it gsave 0.7 setlinewidth pw 0 rlineto 0 th neg rlineto pw neg 0 rlineto closepath stroke grestore % stk: ct rt lt x v get y v get th sub 1 add moveto %%IncludeResource: font Helvetica fHelvetica fnfs 0.8 mul scalefont setfont % 3. The left title gsave dup stringwidth pop fnfs 0.8 mul add exch % leave space took on stack fnfs 0.8 mul hm rmoveto show % left title grestore exch % stk: ct ltw rt % 4. the right title gsave dup stringwidth pop fnfs 0.8 mul add exch % leave space took on stack dup pw exch stringwidth pop fnfs 0.8 mul add sub hm rmoveto show % right title grestore % stk: ct ltw rtw % 5. the center title gsave pw 3 1 roll % stk: ct pw ltw rtw 3 copy % Move to the center of the left room sub add 2 div hm rmoveto % What is the available space in here? add sub fnfs 0.8 mul sub fnfs 0.8 mul sub % stk: ct space_left %%IncludeResource: font Helvetica-Bold fHelvetica-Bold fnfs scalefont setfont cfshow grestore } bind def % Function border: prints virtual page border /border { %def gsave % print four sides 0 setgray x v get y v get moveto 0.7 setlinewidth % of the square pw 0 rlineto 0 ph neg rlineto pw neg 0 rlineto closepath stroke grestore } bind def % Function water: prints a water mark in background /water { %def gsave scx scy moveto rotate %%IncludeResource: font Times-Bold fTimes-Bold 100 scalefont setfont .97 setgray dup stringwidth pop 2 div neg -50 rmoveto show grestore } bind def % Function rhead: prints the right header /rhead { %def lx ly moveto fHelvetica fnfs 0.8 mul scalefont setfont l-show } bind def % Function footer (cf rf lf -> -) /footer { fHelvetica fnfs 0.8 mul scalefont setfont dx dy moveto show snx sny moveto l-show fnx fny moveto c-show } bind def %%EndResource %%BeginResource: procset a2ps-black+white-Prolog 2.0 1 % Function T(ab), jumps to the n-th tabulation in the current line /T { cw mul x0 add y0 moveto } bind def % Function n: move to the next line /n { %def /y0 y0 bfs sub store x0 y0 moveto } bind def % Function N: show and move to the next line /N { Show /y0 y0 bfs sub store x0 y0 moveto } bind def /S { Show } bind def /p { false UL false BX fCourier bfs scalefont setfont Show } bind def /sy { false UL false BX fSymbol bfs scalefont setfont Show } bind def /k { false UL false BX fCourier-Oblique bfs scalefont setfont Show } bind def /K { false UL false BX fCourier-Bold bfs scalefont setfont Show } bind def /c { false UL false BX fCourier-Oblique bfs scalefont setfont Show } bind def /C { false UL false BX fCourier-BoldOblique bfs scalefont setfont Show } bind def /l { false UL false BX fHelvetica bfs scalefont setfont Show } bind def /L { false UL false BX fHelvetica-Bold bfs scalefont setfont Show } bind def /str{ false UL false BX fTimes-Roman bfs scalefont setfont Show } bind def /e{ false UL true BX fHelvetica-Bold bfs scalefont setfont Show } bind def %%EndResource %%EndProlog %%BeginSetup %%IncludeResource: font Courier %%IncludeResource: font Courier-Oblique %%IncludeResource: font Courier-Bold %%IncludeResource: font Times-Roman %%IncludeResource: font Symbol %%IncludeResource: font Courier-BoldOblique %%BeginResource: encoding ISO-8859-1Encoding /ISO-8859-1Encoding [ /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /space /exclam /quotedbl /numbersign /dollar /percent /ampersand /quoteright /parenleft /parenright /asterisk /plus /comma /minus /period /slash /zero /one /two /three /four /five /six /seven /eight /nine /colon /semicolon /less /equal /greater /question /at /A /B /C /D /E /F /G /H /I /J /K /L /M /N /O /P /Q /R /S /T /U /V /W /X /Y /Z /bracketleft /backslash /bracketright /asciicircum /underscore /quoteleft /a /b /c /d /e /f /g /h /i /j /k /l /m /n /o /p /q /r /s /t /u /v /w /x /y /z /braceleft /bar /braceright /asciitilde /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /.notdef /space /exclamdown /cent /sterling /currency /yen /brokenbar /section /dieresis /copyright /ordfeminine /guillemotleft /logicalnot /hyphen /registered /macron /degree /plusminus /twosuperior /threesuperior /acute /mu /paragraph /bullet /cedilla /onesuperior /ordmasculine /guillemotright /onequarter /onehalf /threequarters /questiondown /Agrave /Aacute /Acircumflex /Atilde /Adieresis /Aring /AE /Ccedilla /Egrave /Eacute /Ecircumflex /Edieresis /Igrave /Iacute /Icircumflex /Idieresis /Eth /Ntilde /Ograve /Oacute /Ocircumflex /Otilde /Odieresis /multiply /Oslash /Ugrave /Uacute /Ucircumflex /Udieresis /Yacute /Thorn /germandbls /agrave /aacute /acircumflex /atilde /adieresis /aring /ae /ccedilla /egrave /eacute /ecircumflex /edieresis /igrave /iacute /icircumflex /idieresis /eth /ntilde /ograve /oacute /ocircumflex /otilde /odieresis /divide /oslash /ugrave /uacute /ucircumflex /udieresis /yacute /thorn /ydieresis ] def %%EndResource % Initialize page description variables. /sh 612 def /sw 792 def /llx 24 def /urx 768 def /ury 588 def /lly 24 def /#copies 1 def /th 15.000000 def /fnfs 11 def /bfs 7.493857 def /cw 4.496314 def % Dictionary for ISO-8859-1 support /iso1dict 8 dict begin /fCourier ISO-8859-1Encoding /Courier reencode_font /fCourier-Bold ISO-8859-1Encoding /Courier-Bold reencode_font /fCourier-BoldOblique ISO-8859-1Encoding /Courier-BoldOblique reencode_font /fCourier-Oblique ISO-8859-1Encoding /Courier-Oblique reencode_font /fHelvetica ISO-8859-1Encoding /Helvetica reencode_font /fHelvetica-Bold ISO-8859-1Encoding /Helvetica-Bold reencode_font /fTimes-Bold ISO-8859-1Encoding /Times-Bold reencode_font /fTimes-Roman ISO-8859-1Encoding /Times-Roman reencode_font currentdict end def /bgcolor [ 0 0 0 ] def /bg false def /ul false def /bx false def % The font for line numbering /f# /Helvetica findfont bfs .6 mul scalefont def /fSymbol /Symbol findfont def /hm fnfs 0.25 mul def /pw cw 81.400000 mul def /ph 522.321860 th add def /pmw urx llx sub pw 2 mul sub 1 div def /pmh 0 def /v 0 def /x [ 0 dup pmw add pw add ] def /y [ pmh ph add 0 mul ph add dup ] def /scx sw 2 div def /scy sh 2 div def /snx urx def /sny lly 2 add def /dx llx def /dy sny def /fnx scx def /fny dy def /lx snx def /ly ury fnfs 0.8 mul sub def /sx 0 def /tab 8 def /x0 0 def /y0 0 def %%EndSetup %%Page: (1-2) 1 %%BeginPageSetup /pagesave save def sh 0 translate 90 rotate %%EndPageSetup iso1dict begin gsave llx lly 12 add translate /v 0 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) p n (!) N (! Copyright D.J. Jeffery, 2006jan01.) N (!) N (! plot.f is for making plots of supernova models.) N (!) N (! References:) N (!) N (! \\bibitem[Metcalf et al.\(2004\)]{metcalf2004} Metcalf, M., Reid, J., ) N (! \\&~Cohen, M. 2004,) N (! Fortran 95/2003 Explained) N (! \(Oxford: Oxford University Press\) \(MRC\)) N (! % Pretty good reference book.) N (!) N (! \\bibitem[Nomoto et al.\(1984\)]{nomoto1984}) N (! Nomoto, K., Thielemann, F.-K., \\&~Yokoi, \\apj, 286, 644) N (! % The original W7 paper.) N (! % Here the iron peak yield is .86 solar masses) N (!) N (! \\bibitem[Press et al.\(1992\)]{press1992} Press, W. H., Teukolsky, S. A.,) N (! Vetterling, W. T., \\& Flannery, B. P. 1992,) N (! Numerical Recipes in Fortran) N (! \(Cambridge: Cambridge University Press\)) N (! % Bill should have given me a complimentary copy when) N (! % I co-taught with him in 1993!) N (!) N (! \\bibitem[Thielemann et al.\(1986\)]{thielemann1986} Thielemann, F.-K., ) N (! Nomoto, K.,) N (! \\&~Yokoi, A\\&A, 158, 17) N (! % The final W7 composition paper.) N (! % I sum the iron-peak to be .79 solar masses.) N (! % This is less than the old W7 which had .86 solar masses ) N (! % of iron-peak) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N ( program plot ) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod ! The module is in /aalib/constant.f.) N ( use model_mod ! The module is in /model.) N ( implicit none) N (!) N ( real \(kind=nsngl\) :: abunds\(ncell,natom\)) N ( real \(kind=nsngl\) :: radios\(ncell,nradio_isobar,nradio\)) N (!) N ( integer, parameter :: nfile=10 ) N ( integer, parameter :: nvel_label=50 ) N (!) N ( character :: atmp*80) N ( character :: file_abund*180='abund.ps' ) N ( character :: file_rho*180='rho.ps') N ( integer :: i,j,k,l,m,n ) N ( integer :: iatoms\(nfile\)) N ( integer :: icells\(nfile\)) N ( integer :: ipost\(2\)=1) N ( integer :: iradios\(nfile\)) N ( integer :: iradio_isobars=2) N ( integer :: itmp) N ( real \(kind=nsngl\) :: rhos\(ncell,nfile\)) N ( real \(kind=nsngl\) :: rhocs\(nfile\)) N ( real \(kind=nsngl\) :: sum1,sum2) N ( real \(kind=nsngl\) :: tmp) N ( real \(kind=nprecision\) :: tmpdble) N ( real \(kind=nsngl\) :: tmp1,tmp2,tmp3) N ( real \(kind=nprecision\) :: vefolds\(nfile\)) N ( real \(kind=nsngl\) :: vels\(ncell,nfile\)) N (!) N (tmpa) (Page 1/14) (Jun 12, 07 11:57) title border /v 1 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto ( real \(kind=nsngl\) :: abundmax) p n ( real \(kind=nsngl\) :: abundmin) N ( character :: afile\(nfile\)*180) N ( character :: caption\(nfile\)*80) N ( integer :: iabund=1) N ( integer :: idash\(nfile\) ) N ( integer :: ifile) N ( integer :: ifile_expon\(nfile\) ) N ( integer :: ifile_type\(nfile\) ) N ( integer :: ivel_label=0) N ( integer :: ivel_label_radio=0) N ( integer :: ivel_labelb\(nvel_label\)=0) N ( real \(kind=nsngl\) :: rhomax ) N ( real \(kind=nsngl\) :: rhomin) N ( character :: title*80 ) N ( real \(kind=nsngl\) :: vel_label\(nvel_label\)) N ( real \(kind=nsngl\) :: vel_label_radio\(nvel_label\)) N ( real \(kind=nsngl\) :: velmax) N ( real \(kind=nsngl\) :: velmax_abund) N ( real \(kind=nsngl\) :: velmin) N ( real \(kind=nsngl\) :: velmin_abund) N ( real \(kind=nsngl\) :: xcorr=1.) N ( real \(kind=nsngl\) :: xytable\(2\)=\(/1.e+4,9.e-9/\) ) N ( real \(kind=nsngl\) :: ycorr=1.5 ) N (!) N ( namelist/plotpar/ & ) N ( & abundmax & ) N ( & ,abundmin & ) N ( & ,afile & ) N ( & ,caption & ) N ( & ,file_abund & ) N ( & ,file_rho & ) N ( & ,iabund &) N ( & ,idash & ) N ( & ,ifile & ) N ( & ,ifile_expon & ) N ( & ,ifile_type & ) N ( & ,ipost & ) N ( & ,iradio & ) N ( & ,ivel_center & ) N ( & ,ivel_kms & ) N ( & ,ivel_label & ) N ( & ,ivel_labelb & ) N ( & ,ivel_label_radio & ) N ( & ,rhomax & ) N ( & ,rhomin & ) N ( & ,title & ) N ( & ,vel_label & ) N ( & ,vel_label_radio & ) N ( & ,velmax & ) N ( & ,velmax_abund & ) N ( & ,velmin & ) N ( & ,velmin_abund & ) N ( & ,xytable ) N ( & ) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Before reading the atomic data from /aalib/atomic.') N (!-----------------------------------------------------------------------) N (!) N ( open\(unit=1,file='../../../../aalib/atomic.dat', &) N ( & status='old',action='read'\)) N ( read\(1,atomic_data\)) N ( close\(unit=1\)) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Before reading the plot parameters.') N (tmpa) (Page 2/14) (Jun 12, 07 11:57) title border grestore (Printed by David Jeffery) rhead (tmpa) (1/7) (Tuesday June 12, 2007) footer end % of iso1dict pagesave restore showpage %%Page: (3-4) 2 %%BeginPageSetup /pagesave save def sh 0 translate 90 rotate %%EndPageSetup iso1dict begin gsave llx lly 12 add translate /v 0 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto (!-----------------------------------------------------------------------) p n (!) N ( open\(unit=1,file='plot.par',status='old',action='read'\)) N ( read\(1,plotpar\)) N (! write\(1,plotpar\)) N ( close\(unit=1\)) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Before reading the data files.') N ( write\(*,*\) 'Note we can plot many density profiles, but only') N ( write\(*,*\) 'one abundance file as specified by iabund.') N (!-----------------------------------------------------------------------) N (!) N ( do402 : do i=1,ifile) N ( if\(ifile_type\(i\) .eq. 1\) then) N ( call w7\(trim\(afile\(i\)\)\)) N ( else if\(ifile_type\(i\) .eq. 2\) then) N ( call bravo\(trim\(afile\(i\)\)\)) N ( else if\(ifile_type\(i\) .eq. 3\) then) N ( write\(*,*\) 'Before subroutine hoeflich') N ( call hoeflich\(trim\(afile\(i\)\)\)) N ( end if) N ( write\(*,*\) 'After the model read-in.') N ( write\(*,*\) 'i,ifile,iatom,icell,iradio') N ( write\(*,*\) i,ifile,iatom,icell,iradio) N ( write\(*,*\) 'rhoc,vefold') N ( write\(*,*\) rhoc,vefold) N (!) N ( iatoms\(i\)=iatom) N ( icells\(i\)=icell) N ( iradios\(i\)=iradio) N ( rhos\(:icell,i\)=real\(rho\(:icell\),nsngl\)) N ( rhocs\(i\)=real\(rhoc,nsngl\)) N ( vels\(:icell,i\)=real\(vel\(:icell\),nsngl\)) N ( vefolds\(i\)=real\(vefold,nsngl\)) N (!) N ( rhos\(:icell,i\)=log10\(rhos\(:icell,i\)\)) N ( rhocs\(i\)=log10\(rhocs\(i\)\)) N (!) N ( write\(*,*\) 'After logging density.') N ( if\(i .eq. iabund\) then) N ( abunds\(:icells\(i\),:\)= &) N ( & log10\( real\(max\(abund\(:icell,:\),1.d-30\),nsngl\) \)) N ( radios\(:icells\(i\),:,:\)= &) N ( & log10\( real\(max\(radio\(:icell,:,:\),1.d-30\),nsngl\) \)) N ( end if) N ( end do do402) N (!) N ( rhomax=log10\(rhomax\)) N ( rhomin=log10\(rhomin\)) N ( xytable\(2\)=log10\(xytable\(2\)\)) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Before density plot.') N (!-----------------------------------------------------------------------) N (!) N ( if\(ipost\(1\) .eq. 1\) then) N ( call sm_device\( &) N ( & 'postencap '//file_rho\) ! Square plot) N ( else if\(ipost\(1\) .eq. 2\) then) N ( call sm_device\( &) N ( & 'postportencap '//file_rho\) ! Portrait plot) N ( else) N ( call sm_device\( &) N ( & 'postlandencap '//file_rho\) ! Landscape plot) N ( end if) N ( call sm_graphics) N (tmpa) (Page 3/14) (Jun 12, 07 11:57) title border /v 1 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto ( call sm_defvar\('TeX_strings','1'\)) p n (!) N ( call sm_ticksize\(1.,0.,-1.,10.\)) N (!) N ( print*,'Before sm_limits.') N ( write\(*,*\) 'velmin,velmax,rhomin,rhomax') N ( write\(*,*\) velmin,velmax,rhomin,rhomax ) N ( call sm_limits\(velmin,velmax,rhomin,rhomax\)) N (! print*,'Before sm_box.') N ( call sm_box\(1,2,0,0\)) N (! print*,'Before sm_expand.') N ( call sm_expand\(1.\)) N ( call sm_xlabel\('\\rm Velocity \(km s^{-1}\)'\)) N (! print*,'Before sm_label.') N ( call sm_ylabel\('\\rm Density at 1 day past explosion \(g cm^{-3}\)'\)) N (! call sm_expand\(1.0\)) N (!) N ( print*,'Before sm_grelocate.') N ( call sm_grelocate\(2000,0\) ! Location in screen coordinates: ) N (! ! they run 0--32787) N ( call sm_putlabel\(9,trim\(title\)//'\\rm Density Profile'\)) N (!) N ( write\(*,*\) 'Before connecting the curves.') N ( do410 : do i=1,ifile) N (! write\(*,*\) vels\(:,i\)) N (! write\(*,*\) rhos\(:,i\)) N ( call sm_ltype\(idash\(i\)\)) N ( call sm_conn\(vels\(:,i\),rhos\(:,i\),icells\(i\)\)) N ( if\(ifile_expon\(i\) .eq. 1\) then ! Must update for multiple density pr) N (ofiles. ) N ( itmp=ifile+1) N ( call sm_ltype\(1\)) N ( idash\(itmp\)=1) N ( tmp1=rhocs\(i\)+log10\(exp\(-velmin/vefolds\(i\)\)\)) N ( tmp2=rhocs\(i\)+log10\(exp\(-velmax/vefolds\(i\)\)\)) N ( write\(*,*\) 'tmp1,tmp2') N ( write\(*,*\) tmp1,tmp2) N ( call sm_relocate\(velmin,tmp1\)) N ( call sm_draw\(velmax,tmp2\)) N ( write\(atmp,'\(f5.0\)'\) vefolds\(i\)) N (! caption\(itmp\)=trim\(caption\(2\)\)//trim\(atmp\)//' km s^{-1}') N ( caption\(itmp\)='\\rm Exponential fit, v_{e}=' ) N ( &) N ( & //trim\(atmp\)//' km s^{-1}') N ( end if) N ( end do do410) N ( write\(*,*\) 'After connecting the curves.') N (!) N ( if\(maxval\(ifile_expon\) .eq. 0\) then) N ( itmp=ifile) N ( else) N ( itmp=ifile+1 ! At the moment I allow only one exponential fit.) N ( end if) N ( call sim_scale\(velmin,rhomin,velmax,rhomax\)) N ( call sim_table\(xytable,itmp,idash,caption,xcorr,ycorr\)) N (!) N ( call sm_alpha) N ( call sm_hardcopy) N (!) N ( write\(*,*\) 'The density plot is completed.') N (!) N (!) N ( if\(iabund .ne. 0\) then ! start iabund if) N (!) N ( abundmax=log10\(abundmax\) ) N ( abundmin=log10\(abundmin\) ) N () N (!) N (!-----------------------------------------------------------------------) N (tmpa) (Page 4/14) (Jun 12, 07 11:57) title border grestore (Printed by David Jeffery) rhead (tmpa) (2/7) (Tuesday June 12, 2007) footer end % of iso1dict pagesave restore showpage %%Page: (5-6) 3 %%BeginPageSetup /pagesave save def sh 0 translate 90 rotate %%EndPageSetup iso1dict begin gsave llx lly 12 add translate /v 0 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto ( write\(*,*\)) p n ( write\(*,*\) 'Before abundance plot for file ',iabund,afile\(iabund\)) N ( write\(*,*\) 'At the moment I allow only one model abundances to') N ( write\(*,*\) 'plotted at at time.') N (!-----------------------------------------------------------------------) N (!) N ( if\(ipost\(2\) .eq. 1\) then) N ( call sm_device\( &) N ( & 'postencap '//file_abund\) ! Square plot) N ( else if\(ipost\(2\) .eq. 2\) then) N ( call sm_device\( &) N ( & 'postportencap '//file_abund\) ! Portrait plot) N ( else) N ( call sm_device\( &) N ( & 'postlandencap '//file_abund\) ! Landscape plot) N ( end if) N ( call sm_graphics) N ( call sm_defvar\('TeX_strings','1'\)) N (!) N ( call sm_ticksize\(1.,0.,-1.,10.\)) N (!) N ( print*,'Before sm_limits.') N ( write\(*,*\) 'velmin_abund,velmax_abund,abundmin,abundmax') N ( write\(*,*\) velmin_abund,velmax_abund,abundmin,abundmax) N ( call sm_limits\(velmin_abund,velmax_abund,abundmin,abundmax\)) N (! print*,'Before sm_box.') N ( call sm_box\(1,2,0,0\)) N (! print*,'Before sm_expand.') N ( call sm_expand\(1.\)) N ( call sm_xlabel\('\\rm Velocity \(km s^{-1}\)'\)) N (! print*,'Before sm_label.') N ( call sm_ylabel\('\\rm Mass Fraction'\)) N (! call sm_expand\(1.0\)) N (!) N ( print*,'Before sm_grelocate.') N ( call sm_grelocate\(2000,0\) ! Location in screen coordinates: ) N (! ! they run 0--32787) N ( call sm_putlabel\(9,trim\(title\)//'\\rm Composition Structure'\)) N (!) N ( write\(*,*\) 'Before connecting the curves.') N (! write\(*,*\) iatoms\(iabund\),abunds\(:,26\)) N ( call sm_ltype\(0\)) N ( do420 : do i=1,iatoms\(iabund\)) N ( call sm_conn\(vels\(:,iabund\),abunds\(:,i\),icells\(iabund\)\)) N ( do422 : do j=1,ivel_label) N ( if\(ivel_labelb\(j\) .lt. 0\) then) N ( if\(ivel_labelb\(j\) .eq. -i\) cycle do422 ! Omit one label.) N ( else if\(ivel_labelb\(j\) .gt. 0\) then) N ( if\(ivel_labelb\(j\) .ne. i\) cycle do422 ! Omit all but one.) N ( end if) N ( call interpolation\(icells\(iabund\),) N ( & real\(vels\(:,iabund\),nprecision\), &) N ( & real\(abunds\(:,i\),nprecision\), &) N ( & 1,1,real\(vel_label\(j\),nprecision\),k,tmpdble\)) N ( tmp=real\(tmpdble,nsngl\)) N (! if\(i .eq. 8\) print*,vel_label\(j\),vel\(k,iabund\),abund\(k,i\)) N ( if\(tmp .gt. abundmax .or. tmp .lt. abundmin\) cycle do422) N ( if\(vel_label\(j\) .gt. velmax_abund .or. &) N ( & vel_label\(j\) .lt. velmin_abund\) cycle do422) N ( call sm_relocate\(vel_label\(j\),tmp\)) N ( call sm_putlabel\(5,trim\(aname\(i\)\)\)) N ( end do do422 ) N ( end do do420) N ( call sm_ltype\(1\)) N ( do430 : do i=1,iradios\(iabund\)) N ( do440 : do j=1,iradio_isobars ! Only the first two isotopes in the c) N (hain so far. ) N ( call sm_conn\(vels\(:,iabund\),radios\(:,j,i\),icells\(iabund\)\)) N (!) N (tmpa) (Page 5/14) (Jun 12, 07 11:57) title border /v 1 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto ( do442 : do k=1,ivel_label_radio) p n ( call interpolation\(icells\(iabund\),) N ( & real\(vels\(:,iabund\),nprecision\), &) N ( & real\(radios\(:,j,i\),nprecision\), &) N ( & 1,1,real\(vel_label_radio\(k\),nprecision\),m,tmpdble\)) N ( tmp=real\(tmpdble,nsngl\)) N (! if\(i .eq. 8\) print*,vel_label\(k\),vel\(m,iabund\),radio\(m,j,i\)) N ( if\(tmp .gt. abundmax .or. tmp .lt. abundmin\) cycle do442) N ( if\(vel_label_radio\(k\) .gt. velmax_abund .or. ) N ( &) N ( & vel_label_radio\(k\) .lt. velmin_abund\) cycle do442) N ( call sm_relocate\(vel_label_radio\(k\),tmp\)) N ( call sm_putlabel\(5,trim\(aradio\(j,i\)\)\)) N ( end do do442) N (!) N ( end do do440) N ( end do do430) N ( write\(*,*\) 'After connecting the curves.') N (!) N (! call sim_scale\(velmin_abund,rhomin,velmax_abund,rhomax\)) N (! call sim_table\(xytable,itmp,idash,caption,xcorr,ycorr\)) N (!) N ( call sm_alpha) N ( call sm_hardcopy) N (!) N ( write\(*,*\) 'The abundance plot is completed.') N (!) N ( end if ! end iabund if) N (!) N ( end program plot ) N (!) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (! Copyright D.J. Jeffery, 2006jan01.) N (!) N (! This is a set of subroutines for reading in supernova models.) N (! They come in all different formats.) N (!) N (! References:) N (!) N (! \\bibitem[Metcalf et al.\(2004\)]{metcalf2004} Metcalf, M., Reid, J., \\&~Cohe) N (n, M. 2004,) N (! Fortran 95/2003 Explained) N (! \(Oxford: Oxford University Press\) \(MRC\)) N (! % Pretty good reference book.) N (!) N (! \\bibitem[Nomoto et al.\(1984\)]{nomoto1984}) N (! Nomoto, K., Thielemann, F.-K., \\&~Yokoi, \\apj, 286, 644) N (! % The original W7 paper.) N (! % Here the iron peak yield is .86 solar masses) N (!) N (! \\bibitem[Press et al.\(1992\)]{press1992} Press, W. H., Teukolsky, S. A.,) N (! Vetterling, W. T., \\& Flannery, B. P. 1992,) N (! Numerical Recipes in Fortran) N (! \(Cambridge: Cambridge University Press\)) N (! % Bill should have given me a complimentary copy when) N (! % I co-taught with him in 1993!) N (!) N (! \\bibitem[Thielemann et al.\(1986\)]{thielemann1986} Thielemann, F.-K., Nomo) N (to, K.,) N (! \\&~Yokoi, A\\&A, 158, 17) N (! % The final W7 composition paper.) N (! % I sum the iron-peak to be .79 solar masses.) N (! % This is less than the old W7 which had .86 solar masses of iron-) N (peak) N (!) N (tmpa) (Page 6/14) (Jun 12, 07 11:57) title border grestore (Printed by David Jeffery) rhead (tmpa) (3/7) (Tuesday June 12, 2007) footer end % of iso1dict pagesave restore showpage %%Page: (7-8) 4 %%BeginPageSetup /pagesave save def sh 0 translate 90 rotate %%EndPageSetup iso1dict begin gsave llx lly 12 add translate /v 0 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto (! \\bibitem[Yoon \\&~Langer\(2005\)]{yoon2005}) p n (! Yoon, S.-C., \\&~Langer, N. 2005, A\\&A, 435, 985) N (! % On the evolution of rapidly rotating massive white dwarfs toward supe) N (rnovae) N (! % or collapses.) N (! % On p. 977 they give a figure of binding energy for massive rotating W) N (Ds.) N (! % They also give fitting formula for those binding energies.) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N ( module model_mod) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod ! The module is in /aalib/constant.f) N ( implicit none) N (!) N ( integer, parameter :: ncell=400) N ( real \(kind=nprecision\) :: abund\(ncell,natom\)=0.d0) N ( integer :: iatom ! Model upper limit on atomic number.) N ( integer :: icell ! Number of cells in model.) N ( integer :: iradio) N ( integer :: ivel_center=0 ! 0 for none \(default\); 1 for simple averaging;) N (! ! 2 for tricky, but it's about the same as 1.) N ( integer :: ivel_kms=0 ! 0 for no conversion; 1 for conversion to km/s) N (.) N ( real \(kind=nprecision\) :: radio\(ncell,nradio_isobar,nradio\)=0.d0) N ( real \(kind=nprecision\) :: rho\(ncell\)) N ( real \(kind=nprecision\) :: rhoc) N ( real \(kind=nprecision\) :: time ! Time of model in seconds. ) N ( real \(kind=nprecision\) :: vefold) N ( real \(kind=nprecision\) :: vel\(ncell\)) N ( real \(kind=nprecision\) :: xke_erg ) N ( real \(kind=nprecision\) :: xmass_cell\(ncell\) ) N ( real \(kind=nprecision\) :: xmass_cell_del\(ncell\) ) N ( real \(kind=nprecision\) :: xmass_gram) N (!) N ( end module model_mod) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (!) N (!) N ( subroutine model_diagnostic) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod ! The module is in /aalib/constant.f) N ( use model_mod) N ( implicit none) N (!) N ( integer :: i,j,k,l,m,n ) N ( real \(kind=nprecision\) :: vel_ipe) N ( real \(kind=nprecision\) :: xmass_ipe) N ( real \(kind=nprecision\) :: xmass_nico ) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Getting the exponential fitted vefold and rhoc.') N (!-----------------------------------------------------------------------) N (!) N ( vefold=sqrt\(sixth*xke_erg/xmass_gram\) ! Exponential model fit for vef) N (old.) N ( rhoc=xmass_gram/\(8.d0*pi*\(vefold*daysec\)**3\)) N (! ! Exponential model fit for rhoc.) N ( write\(*,*\) 'Fitted vefold and rhoc for 1 day.') N ( write\(*,*\) 'vefold,rhoc,xmass_gram,xke_erg') N ( write\(*,'\(1p,4e25.15\)'\) vefold,rhoc,xmass_gram,xke_erg) N (tmpa) (Page 7/14) (Jun 12, 07 11:57) title border /v 1 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto (!) p n (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Checking the Ni/Co-56, IPE mass, and v_ipe.') N (!-----------------------------------------------------------------------) N (!) N (! Sum up the Ni/Co-56 mass.) N ( xmass_nico=sum\( \( radio\(:icell,1,1\)+radio\(:icell,2,1\) \) &) N ( & *xmass_cell_del\(:icell\) \)) N (!) N (! Sum up the IPE mass.) N ( xmass_ipe=0.d0) N ( do410 : do i=1,icell) N ( xmass_ipe=xmass_ipe+sum\(abund\(i,21:32\)\)*xmass_cell_del\(i\)) N ( end do do410) N (!) N (! Determine the fiducial IPE mass core.) N (! Usually the velocities and masses are both at the cell walls I think.) N ( call interpolation\(icell,xmass_cell\(:icell\), &) N ( & vel\(:icell\),1,1,xmass_ipe,k,vel_ipe\)) N (!) N ( write\(*,*\) 'Ni/Co56 mass,IPE mass,g-factor,v_ipe,vel\(k\)') N ( write\(*,*\) xmass_nico,xmass_ipe,xmass_nico/xmass_ipe, &) N ( & vel_ipe,vel\(k\)) N ( write\(*,*\) 'Si,S,IPE,Ni-56 at vel\(k\)') N ( write\(*,*\) abund\(k,14\),abund\(k,16\),sum\(abund\(k,21:28\)\), &) N ( & radio\(k,1,1\),vel\(k\)) N ( write\(*,*\) 'Si,S,IPE at vel\(k+1\)') N ( write\(*,*\) abund\(k+1,14\),abund\(k+1,16\),sum\(abund\(k+1,21:28\)\),) N ( & radio\(k+1,1,1\),vel\(k+1\)) N (!) N ( end subroutine model_diagnostic) N (!) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (! Subroutine bravo is for reading the model files of Bravo & Garcia-Senz. ) N (!) N ( subroutine bravo\(afile\)) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod) N ( use model_mod) N ( implicit none) N (!) N ( integer, parameter :: nbravo=15) N (!) N ( character, intent\(in\) :: afile*\(*\)) N (!) N ( integer :: i,j,k,l,m,n) N ( integer :: ibravo\(nbravo\)=\(/ 2, 6, 8,10, 12,14,16,18, &) N ( & 20,24,25,26, 27,28,30/\) ) N ( real \(kind=nprecision\) :: tmp) N ( real \(kind=nprecision\) :: tmp1,tmp2,tmp3) N (!) N ( namelist/param/iatom,time) N (!) N ( open\(unit=1,file='/home/jeffery/public_html/aalib/atomic.dat', &) N ( & status='old',action='read'\)) N ( read\(1,atomic_data\)) N ( close\(unit=1\)) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'bravo: Before zeroing and reading in.') N (!-----------------------------------------------------------------------) N (!) N ( abund\(:,:\)=0.d0) N (tmpa) (Page 8/14) (Jun 12, 07 11:57) title border grestore (Printed by David Jeffery) rhead (tmpa) (4/7) (Tuesday June 12, 2007) footer end % of iso1dict pagesave restore showpage %%Page: (9-10) 5 %%BeginPageSetup /pagesave save def sh 0 translate 90 rotate %%EndPageSetup iso1dict begin gsave llx lly 12 add translate /v 0 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto ( radio\(:,:,:\)=0.d0) p n ( xke_erg=0.d0) N ( xmass_cell_del\(:\)=0.d0) N (!) N ( open\(unit=1,file=afile,status='old',action='read'\) ! Open bravo data file.) N ( read\(1,param\)) N ( call readcopy\(1,0,0,' ',icell\)) N ( write\(*,*\)) N ( write\(*,"\('icell',3x,'vel \(cm/s\)',4x,'rho \(cgs\)', &) N ( & 2x,'m_i \(m_sun\)',5x,'ratio', &) N ( & 7x,'sum',3x,'X Ni-56'\)"\)) N ( do410 : do i=1,icell) N ( read\(1,*\) j,tmp,xmass_cell\(i\),vel\(i\), &) N ( & \(abund\(i,ibravo\(j\)\),j=1,nbravo\), &) N ( & \(radio\(i,1,j\),radio\(i,2,j\),radio\(i,3,j\),j=1,2\), &) N ( & radio\(i,1,4\),radio\(i,2,4\),radio\(i,3,4\), &) N ( & radio\(i,2,3\),radio\(1,3,3\) ! The 1st nuclide Co-55 is not) N ( in Bravo.) N ( abund\(i,28\)=abund\(i,28\)+radio\(i,1,1\)+radio\(i,1,2\)) N ( abund\(i,27\)=abund\(i,27\)+radio\(i,2,1\)+radio\(i,2,2\)) N ( abund\(i,26\)=abund\(i,26\)+radio\(i,3,1\)+radio\(i,3,2\)) N (!) N ( abund\(i,26\)=abund\(i,26\)+radio\(i,2,3\) ) N ( abund\(i,25\)=abund\(i,25\)+radio\(i,3,3\)) N (!) N ( abund\(i,22\)=abund\(i,22\)+radio\(i,1,4\)) N ( abund\(i,21\)=abund\(i,21\)+radio\(i,2,4\)) N ( abund\(i,20\)=abund\(i,20\)+radio\(i,3,4\)) N (!) N (! Note in radio\(i,j,k\) that the last index is the decay chain, the 2nd to) N (! last is the radionuclide in the decay chain.) N (! Bravo doesn't have Co-55 abundances in his table.) N (!) N ( if\(i .gt. 1\) then) N ( xmass_cell_del\(i\)=xmass_cell\(i\)-xmass_cell\(i-1\)) N ( tmp1=vel\(i-1\)) N ( tmp2=vel\(i-1\)**3) N ( else) N ( xmass_cell_del\(i\)=xmass_cell\(i\)) N ( tmp1=0.d0 ) N ( tmp2=0.d0) N ( end if) N ( rho\(i\)=\(xmass_cell_del\(i\)*xmass_sun\) &) N ( & /\(pi4L3*\(vel\(i\)**3-tmp2\)*daysec**3\)) N ( xke_erg=xke_erg+.5d0*\( .5d0*\(vel\(i\)+tmp1\) \)**2 &) N ( & *xmass_cell_del\(i\)*xmass_sun ) N ( tmp3=sum\(abund\(i,:iatom\)\)) N ( write\(*,910\) i,vel\(i\),rho\(i\), &) N ( & xmass_cell\(i\),tmp/\(vel\(i\)*time\),tmp3, &) N ( & abund\(i,28\)) N ( 910 format\(i5,1p,3e13.4,0p,3f10.4\)) N ( end do do410) N ( close\(unit=1\) ) N (!) N ( xmass_gram=xmass_cell\(icell\)*xmass_sun) N (!) N ( call model_diagnostic) N ( call velocity_centering) N ( call velocity_cms_to_kms) N (!) N ( i=1) N ( write\(*,*\) 'i,iatom,abund\(i,8\),awt\(8\)') N ( write\(*,*\) i,iatom,abund\(i,8\),awt\(8\)) N (!) N ( end subroutine bravo) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (! Subroutine hoeflich is for reading the old model files of Peter Hoeflich. ) N (tmpa) (Page 9/14) (Jun 12, 07 11:57) title border /v 1 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto (!) p n ( subroutine hoeflich\(afile\)) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod) N ( use model_mod) N ( implicit none) N (!) N ( integer, parameter :: nnatom=30) N (!) N ( character, intent\(in\) :: afile*\(*\)) N (!) N ( integer :: i,j,k,l,m,n) N ( integer :: iiatom\(nnatom\)) N ( character :: record*180 ) N ( real \(kind=nprecision\) :: tmp) N ( real \(kind=nprecision\) :: tmp1,tmp2,tmp3) N ( real \(kind=nprecision\) :: xmass_mult ) N (!) N ( namelist/param/iatom,iiatom,icell,rhoc,time,xmass_mult) N (!) N ( open\(unit=1,file='/home/jeffery/public_html/aalib/atomic.dat', &) N ( & status='old',action='read'\)) N ( read\(1,atomic_data\)) N ( close\(unit=1\)) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'hoeflich: Before zeroing and reading in.') N (!-----------------------------------------------------------------------) N (!) N ( abund\(:,:\)=0.d0) N ( radio\(:,:,:\)=0.d0) N ( xke_erg=0.d0) N ( xmass_cell_del\(:\)=0.d0) N (!) N ( open\(unit=1,file=afile,status='old',action='read'\) ! Open bravo data file.) N ( read\(1,param\)) N ( call readcopy\(1,0,0,' ',k\)) N ( read\(1,*\)) N ( read\(1,*\) record) N ( write\(*,*\) trim\(record\)) N ( write\(*,"\('icell',3x,'vel \(cm/s\)',4x,'rho \(cgs\)', &) N ( & 2x,'m_i \(m_sun\)',5x,'ratio', &) N ( & 7x,'sum',3x,'X Ni-56'\)"\)) N ( do410 : do i=1,icell) N ( read\(1,*\) j,xmass_cell\(i\),rho\(i\),vel\(i\),tmp2,tmp3,radio\(i,1,1\)) N ( end do do410) N ( do420 : do i=1,11) N ( read\(1,*\)) N ( end do do420) N ( do430 : do i=1,icell) N ( read\(1,*\) j) N ( write\(*,*\) i,j) N ( read\(1,*\) \(abund\(i,iiatom\(j\)\),j=1,7\)) N ( read\(1,*\) \(abund\(i,iiatom\(j\)\),j=8,14\)) N ( read\(1,*\) \(abund\(i,iiatom\(j\)\),j=15,21\)) N ( abund\(i,28\)=abund\(i,28\)+radio\(i,1,1\)) N ( abund\(i,26\)=abund\(i,26\)-radio\(i,1,1\)) N ( write\(*,*\) i,abund\(i,26\),abund\(i,28\),radio\(i,1,1\)) N ( end do do430) N (!) N (! Note in radio\(i,j,k\) that the last index is the decay chain, the 2nd to) N (! last is the radionuclide in the decay chain.) N (!) N ( rhoc=rhoc*\(time/daysec\)**3) N ( rho\(:\)=rho\(:\)*rhoc) N ( xmass_cell\(:\)=xmass_cell\(:\)*xmass_mult) N (tmpa) (Page 10/14) (Jun 12, 07 11:57) title border grestore (Printed by David Jeffery) rhead (tmpa) (5/7) (Tuesday June 12, 2007) footer end % of iso1dict pagesave restore showpage %%Page: (11-12) 6 %%BeginPageSetup /pagesave save def sh 0 translate 90 rotate %%EndPageSetup iso1dict begin gsave llx lly 12 add translate /v 0 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto (!) p n ( do440 : do i=1,icell) N ( if\(i .gt. 1\) then) N ( xmass_cell_del\(i\)=xmass_cell\(i\)-xmass_cell\(i-1\)) N ( tmp1=vel\(i-1\)) N ( tmp2=vel\(i-1\)**3) N ( else) N ( xmass_cell_del\(i\)=xmass_cell\(i\)) N ( tmp1=0.d0) N ( tmp2=0.d0) N ( end if) N ( tmp=\(xmass_cell_del\(i\)*xmass_sun\) &) N ( & /\(pi4L3*\(vel\(i\)**3-tmp2\)*daysec**3\)) N ( xke_erg=xke_erg+.5d0*\( .5d0*\(vel\(i\)+tmp1\) \)**2 &) N ( & *xmass_cell_del\(i\)*xmass_sun) N ( tmp3=sum\(abund\(i,:iatom\)\)) N ( write\(*,910\) i,vel\(i\),rho\(i\), &) N ( & xmass_cell\(i\),tmp/rho\(i\),tmp3, &) N ( & abund\(i,28\)) N ( 910 format\(i5,1p,3e13.4,0p,3f10.4\)) N ( end do do440) N ( close\(unit=1\)) N (!) N ( xmass_gram=xmass_cell\(icell\)*xmass_sun) N (!) N ( call model_diagnostic) N ( call velocity_centering) N ( call velocity_cms_to_kms) N (!) N ( i=1) N ( write\(*,*\) 'i,iatom,abund\(i,8\),awt\(8\)') N ( write\(*,*\) i,iatom,abund\(i,8\),awt\(8\)) N (!) N ( end subroutine hoeflich ) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (! Subroutine w7 is for reading the w7.dat file. Not finished.) N (!) N ( subroutine w7\(afile\)) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod ) N ( use model_mod ) N ( implicit none) N (!) N ( character, intent\(in\) :: afile*\(*\)) N ( real \(kind=nprecision\) :: amass\(nnuclide\) ) N ( real \(kind=nprecision\) :: frac\(nnuclide\)) N ( integer :: i,j,k,l,m,n ) N ( integer :: inuclide\(nnuclide\) ) N ( character :: record*80) N ( real \(kind=nprecision\) :: rho_cell ) N ( real \(kind=nprecision\) :: sum1,sum2,sum3) N ( real \(kind=nprecision\) :: tmp ) N ( real \(kind=nprecision\) :: tmp1,tmp2,tmp3) N (!) N ( real \(kind=nprecision\) :: acarbon=.4875) N ( real \(kind=nprecision\) :: aoxygen=.4875) N ( real \(kind=nprecision\) :: aneon=.025) N ( integer :: icell_comp) N ( integer :: iinuclide) N ( integer :: itmp) N ( real \(kind=nprecision\) :: xmass ) N ( namelist/param/iatom,icell,icell_comp,iinuclide, &) N ( & iradio,time,xmass) N (!) N ( character :: anuclide\(nnuclide\)*4 ) N (tmpa) (Page 11/14) (Jun 12, 07 11:57) title border /v 1 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto ( real \(kind=nprecision\) :: radius\(ncell\)) p n (!) N ( namelist/runs/anuclide,radius,rho,vel,xmass_cell) N (!) N ( open\(unit=1,file='/home/jeffery/public_html/aalib/atomic.dat', &) N ( & status='old',action='read'\)) N ( read\(1,atomic_data\) ) N ( close\(unit=1\)) N (!) N ( open\(unit=1,file=afile,status='old',action='read'\) ! Open w7 data file.) N ( read\(1,param\)) N ( read\(1,runs\)) N (!) N (! Dealing with the first four nuclides.) N ( inuclide\(1:4\)=1 ! All hydrogen or they will be.) N ( amass\(1:2\)=1.d0) N ( amass\(3\)=2.d0) N ( amass\(4\)=3.d0) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'w7: Before identifiying the nuclides with '// &) N ( & 'their atoms.') N (!-----------------------------------------------------------------------) N (!) N ( write\(*,"\('nuclide',' Atomic No.',' Atomic mass'\)"\)) N ( do402 : do i=5,nnuclide ! The first four a N,P,D,T: dealt with above.) N ( do404 : do j=1,natom) N ( if\(anuclide\(i\)\(1:2\) .eq. aname2\(j\)\) then) N ( inuclide\(i\)=j) N ( read\(anuclide\(i\)\(3:4\),'\(f3.0\)'\) amass\(i\) ) N ( end if ) N ( end do do404 ) N ( write\(*,'\(i7,i11,f11.0\)'\) i,inuclide\(i\),amass\(i\)) N ( end do do402) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Before converting density and masses.') N (!-----------------------------------------------------------------------) N (!) N (! Convert to density at day 1.) N ( rho\(:\)=10.d0**rho\(:\)*\(time/daysec\)**3 ) N (! Conversion to solar masses. Masses must be in log of fractional mass.) N ( xmass_cell\(:icell\)=\( 10.d0**xmass_cell\(:icell\) \)*xmass ) N ( xmass_gram=0.d0 ! Mass from the mass in the cells.) N ( xke_erg=0.d0 ! Kinetic energy.) N ( abund\(:,:\)=0.d0) N ( frac\(:\)=0.d0) N ( radio\(:,:,:\)=0.d0) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Before calculating mass and kinetic energy and') N ( write\(*,*\) 'doing some consistency checks which passed so-so.') N ( write\(*,*\)) N ( write\(*,"\(' zone',16x,'rho given',11x,'rho calculated', &) N ( & 1x,'rho ratio',1x,'cell mass'\)"\)) N (! ) N (!-----------------------------------------------------------------------) N (!) N ( do410 : do i=1,icell) N ( if\(i .eq. 1\) then) N ( tmp=0.d0) N ( tmp1=0.d0) N ( tmp2=0.d0) N ( xmass_cell_del\(i\)=xmass_cell\(i\)) N ( else) N ( tmp=vel\(i-1\)**3) N (tmpa) (Page 12/14) (Jun 12, 07 11:57) title border grestore (Printed by David Jeffery) rhead (tmpa) (6/7) (Tuesday June 12, 2007) footer end % of iso1dict pagesave restore showpage %%Page: (13-14) 7 %%BeginPageSetup /pagesave save def sh 0 translate 90 rotate %%EndPageSetup iso1dict begin gsave llx lly 12 add translate /v 0 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto ( tmp1=xmass_cell\(i-1\)) p n ( tmp2=vel\(i-1\)) N ( xmass_cell_del\(i\)=xmass_cell\(i\)-xmass_cell\(i-1\)) N ( end if) N (!) N ( rho_cell=xmass_cell_del\(i\)*xmass_sun &) N ( & /\(pi4L3*\(vel\(i\)**3-tmp\)*daysec**3\)) N ( write\(*,'\(i5,1p,2e25.15,0p,2f10.5\)'\) i,rho\(i\),rho_cell, &) N ( & \(rho\(i\)/rho_cell\), &) N ( & xmass_cell\(i\)) N ( xmass_gram=xmass_gram+rho\(i\)*\(pi4L3*\(vel\(i\)**3-tmp\)*daysec**3\)) N (!) N ( xke_erg=xke_erg+.5d0*\( .5d0*\(vel\(i\)+tmp2\) \)**2 ) N ( &) N ( & *rho\(i\)*\(pi4L3*\(vel\(i\)**3-tmp\)*daysec**3\)) N (!) N ( if\(i .le. icell_comp\) then ! w7 abundances only given up to icell_comp) N ( 110 continue) N ( read\(1,'\(a80\)'\) record) N (! write\(*,*\) record) N ( if\(record\(2:4\) .ne. 'ZO='\) go to 110) N ( do412 : do j=5,iinuclide+1,5 ! Reading the nuclide fractions.) N ( read\(1,'\(2x,5\(4x,e10.3\)\)'\) &) N ( & \(frac\(k\),k=j-4,min\(j,iinuclide\)\) ! reading the nuclide fractions) N ( end do do412 ) N ( sum3=sum\(frac\(:\)*amass\(:\)\) ) N ( frac\(:\)=frac\(:\)*amass\(:\)/sum3) N ( do414 : do j=1,nnuclide) N ( abund\(i,inuclide\(j\)\)=abund\(i,inuclide\(j\)\)+frac\(j\)) N ( end do do414 ) N ( radio\(i,1,1\)=frac\(204\) ! Ni-56 parent.) N ( radio\(i,2,1\)=frac\(195\) ! Co-56 daughter.) N ( else ! Set outer layer abundances.) N ( abund\(i,6\)=acarbon) N ( abund\(i,8\)=aoxygen) N ( abund\(i,10\)=aneon) N ( radio\(i,1,1\)=0.d0) N ( radio\(i,2,1\)=0.d0) N ( end if) N (!) N ( end do do410) N (!) N ( close\(unit=1\) ! Close w7 data file.) N (!) N (!-----------------------------------------------------------------------) N ( write\(*,*\)) N ( write\(*,*\) 'Checking the mass consistency and time') N (!-----------------------------------------------------------------------) N (!) N ( sum1=xmass_gram/xmass_sun) N ( tmp=\(xmass/sum1\)**third*time) N ( write\(*,*\) 'xmass,xmass cal.,time,time cal.,kinetic energy') N ( write\(*,'\(2f25.15,1p,3e25.15\)'\) xmass,sum1,time,tmp,xke_erg) N (!) N ( xmass_gram=xmass*xmass_sun ! This is just for w7 accuracy.) N (!) N ( call model_diagnostic) N ( call velocity_centering) N ( call velocity_cms_to_kms) N (!) N ( end subroutine w7) N (!) N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (! Conversion to velocity centered in mass cell.) N (!) N (! The simple average agrees with the more elaborate result to) N (! better than 1 % almost always for any reasonably well resolved model) N (! based on one example test. But it seems that this should be generally) N (tmpa) (Page 13/14) (Jun 12, 07 11:57) title border /v 1 store /x0 x v get 3.147420 add sx cw mul add store /y0 y v get bfs th add sub store x0 y0 moveto (! true. One can always test if one is unsure.) p n (!) N ( subroutine velocity_centering) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod) N ( use model_mod) N ( implicit none) N (!) N ( integer :: i,j,k,l,m,n) N ( real \(kind=nprecision\) :: tmp,tmp1) N (!) N ( if\(ivel_center .eq. 0\) return) N (!) N ( tmp=0.d0) N ( do410 : do i=1,icell) N ( tmp1=vel\(i\)) N ( if\(ivel_center .eq. 1\) then) N ( vel\(i\)=.5d0*\(tmp+tmp1\)) N (! print*,'vel\(i\)') N (! print*,vel\(i\)) N ( else) N ( vel\(i\)=sqrt\(xmass_cell_del\(i\)*xmass_sun &) N ( & /\(pi4*rho\(i\)*daysec**3*\(tmp1-tmp\)\)\)) N (! print*,vel\(i\)) N ( end if) N ( tmp=tmp1) N ( end do do410) N (!) N ( end subroutine velocity_centering) N () N (!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12) N (!) N (! Conversion of velocity from cm/s to km/s.) N (!) N ( subroutine velocity_cms_to_kms) N ( use numerical) N ( use physical) N ( use astronomical) N ( use atomic_mod) N ( use model_mod) N ( implicit none) N (!) N ( if\(ivel_kms .eq. 0\) return) N (!) N ( vel\(:\)=vel\(:\)*1.d-5 ) N ( vefold=vefold*1.d-5 ) N (!) N ( end subroutine velocity_cms_to_kms) N (tmpa) (Page 14/14) (Jun 12, 07 11:57) title border grestore (Printed by David Jeffery) rhead (tmpa) (7/7) (Tuesday June 12, 2007) footer end % of iso1dict pagesave restore showpage %%Trailer end %%EOF