added original Fortran sources for PROGRAM1.

Author: alsam

examples/scalar_law/PROGRAM1/fortran/GNUmakefile

@@ -0,0 +1,10 @@
+# to run this program after it has been made, type
+#   flinearad > output
+# to see the results after running, type
+#   xmgrace output
+FC=gfortran
+flinearad : linaddmain.o consdiff.o linearad.o riemprob.o upwind.o
+	$(FC) -o $@ linaddmain.o consdiff.o linearad.o riemprob.o \
+	  upwind.o
+clean :
+	rm -f *.o flinearad

examples/scalar_law/PROGRAM1/fortran/consdiff.f

@@ -0,0 +1,66 @@
+c***********************************************************************
+c  Copyright 2006 John A. Trangenstein
+c
+c  This software is made available for research and instructional use 
+c  only. 
+c  You may copy and use this software without charge for these 
+c  non-commercial purposes, provided that the copyright notice and 
+c  associated text is reproduced on all copies.  
+c  For all other uses (including distribution of modified versions), 
+c  please contact the author at
+c    John A. Trangenstein
+c    Department of Mathematics
+c    Duke University
+c    Durham, NC 27708-0320
+c    USA
+c  or
+c    [email protected]
+c  
+c  This software is made available "as is" without any assurance that it
+c  is completely correct, or that it will work for your purposes.  
+c  Use the software at your own risk.
+c***********************************************************************
+      subroutine consdiff(fic,lac,fif,laf,fix,lax,ifirst,ilast,
+     &  x,flux,
+     &  conservd)
+      integer fic,lac,fif,laf,fix,lax,ifirst,ilast
+      double precision
+     &  x(fix:lax+1),
+     &  flux(fif:laf+1)
+      double precision
+     &  conservd(fic:lac)
+      integer ic
+c     ******************************************************************
+c     update conservd to new time
+c     ******************************************************************
+      do ic=ifirst,ilast
+        conservd(ic) = conservd(ic)
+     &               - (flux(ic+1)-flux(ic)) / (x(ic+1)-x(ic))
+      enddo
+      return
+      end
+
+      function stabledt(fc,lc,fm,lm,ifirst,ilast,
+     &  conserved,lambda_cell,x)
+c     ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      include "./const.i"
+c     ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      double precision stabledt
+      integer fc,lc,fm,lm,ifirst,ilast
+      double precision conserved(fc:lc),lambda_cell(fc:lc),x(fm:lm+1)
+      integer ic
+      double precision abs_lambda,dx
+c     ******************************************************************
+c     compute stable timestep
+c     ******************************************************************
+      stabledt=huge
+      do ic=ifirst,ilast
+        abs_lambda=abs(lambda_cell(ic))
+        dx = x(ic+1)-x(ic)
+        if (abs_lambda.gt.roundoff*dx) then
+          stabledt=min(stabledt,dx/abs_lambda)
+        endif
+      enddo
+c     if (stabledt.le.0.) call abort()
+      return
+      end 

examples/scalar_law/PROGRAM1/fortran/const.i

@@ -0,0 +1,57 @@
+c***********************************************************************
+c  Copyright 2006 John A. Trangenstein
+c
+c  This software is made available for research and instructional use 
+c  only. 
+c  You may copy and use this software without charge for these 
+c  non-commercial purposes, provided that the copyright notice and 
+c  associated text is reproduced on all copies.  
+c  For all other uses (including distribution of modified versions), 
+c  please contact the author at
+c    John A. Trangenstein
+c    Department of Mathematics
+c    Duke University
+c    Durham, NC 27708-0320
+c    USA
+c  or
+c    [email protected]
+c  
+c  This software is made available "as is" without any assurance that it
+c  is completely correct, or that it will work for your purposes.  
+c  Use the software at your own risk.
+c***********************************************************************
+c static char sccsid[] = "%W%\t%G%"
+      double precision 
+     &  roundoff,small,huge,undefind,lnrndoff,lnsmall
+      integer ihuge,inhuge
+      common/machine/roundoff,small,huge,undefind,lnrndoff,lnsmall,
+     &  ihuge,inhuge
+      double precision 
+     &  zero,sixtyfourth,thirtysecond,sixteenth,tenth,eighth,sixth,
+     &  fourth,third,point4,half,twothird,pt75,onemert2,rt75,one,
+     &  fourthir,onept5,two,three,pi,four,seven,nine,ten
+      parameter (zero=0.d0)
+      parameter (sixtyfourth=0.015625d0)
+      parameter (thirtysecond=0.03125d0)
+      parameter (sixteenth=0.0625d0)
+      parameter (tenth=0.1d0)
+      parameter (eighth=0.125d0)
+      parameter (sixth=0.16666666666667d0)
+      parameter (fourth=.25d0)
+      parameter (third=.333333333333333d0)
+      parameter (point4=.4d0)
+      parameter (half=.5d0)
+      parameter (twothird=.66666666666667d0)
+      parameter (pt75=.75d0)
+      parameter (onemert2=.75688326556578578920d0)
+      parameter (rt75=.8660254037844d0)
+      parameter (one=1.d0)
+      parameter (fourthir=1.33333333333333d0)
+      parameter (onept5=1.5d0)
+      parameter (two=2.d0)
+      parameter (three=3.d0)
+      parameter (pi=3.14159265358979323846d0)
+      parameter (four=4.d0)
+      parameter (seven=7.d0)
+      parameter (nine=9.d0)
+      parameter (ten=10.d0)

examples/scalar_law/PROGRAM1/fortran/linaddmain.f

@@ -0,0 +1,108 @@
+c***********************************************************************
+c  Copyright 2006 John A. Trangenstein
+c
+c  This software is made available for research and instructional use 
+c  only. 
+c  You may copy and use this software without charge for these 
+c  non-commercial purposes, provided that the copyright notice and 
+c  associated text is reproduced on all copies.  
+c  For all other uses (including distribution of modified versions), 
+c  please contact the author at
+c    John A. Trangenstein
+c    Department of Mathematics
+c    Duke University
+c    Durham, NC 27708-0320
+c    USA
+c  or
+c    [email protected]
+c  
+c  This software is made available "as is" without any assurance that it
+c  is completely correct, or that it will work for your purposes.  
+c  Use the software at your own risk.
+c***********************************************************************
+      program linaddmain
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      integer ncells
+!      parameter (ncells=100)
+      parameter (ncells=10000)
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      include "linearad.i"
+      include "const.i"
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      double precision stabledt
+      external stabledt
+
+      double precision
+     &  u(-2:ncells+1),
+     &  x(0:ncells),
+     &  flux(0:ncells),
+     &  dfdu(-2:ncells+1)
+c    &  flux_cell(-2:ncells+1),
+c    &  uside(0:ncells)
+      integer nsteps
+      double precision tmax,cfl
+
+      integer i,fc,lc,fm,lm,fs,ls,ifirst,ilast
+      integer istep
+      double precision dt,t
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      roundoff=1.d-14
+      small=1.d-20
+      huge=1.d300
+      undefind=1.d300
+      lnrndoff=14.d0
+      lnsmall=20.d0
+
+c     problem-specific parameters:
+      jump=0
+      x_left=-0.2d0
+      x_right=1.d0
+      statelft=2.d0
+      statergt=0.d0
+      velocity=1.d0
+
+      nsteps=10000
+      tmax=0.8d0
+      cfl=0.9d0
+
+c     array bounds for subroutine calls:
+      fc=-2
+      lc=ncells+1
+      fm=0
+      lm=ncells-1
+      fs=0
+      ls=ncells-1
+      ifirst=0
+      ilast=ncells-1
+
+c     initialization:
+      call initsl(ncells,fc,lc,fm,lm,ifirst,ilast, u,x)
+c     do i=0,ncells
+c       print *,"x(",i,") = ",x(i)
+c     enddo
+c     do i=-2,ncells+1
+c       print *,"u(",i,") = ",u(i)
+c     enddo
+      call bcmesh(fm,lm,ncells, x)
+      call bccells(fc,lc,ncells, u)
+      call fluxderv(fc,lc,fc,lc, u, dfdu)
+
+      istep=0
+      t=0.d0
+      do while (istep.lt.nsteps .and. t.lt.tmax)
+        call bccells(fc,lc,ncells, u)
+        call fluxderv(fc,lc,fc,lc, u, dfdu)
+        dt=cfl*stabledt(fc,lc,fm,lm,ifirst,ilast, u, dfdu,x)
+        call method(dt,fc,lc,fm,lm,fs,ls,ifirst,ilast, u, flux)
+        call consdiff(fc,lc,fs,ls,fm,lm,ifirst,ilast, x,flux, u)
+        t=t+dt
+        istep=istep+1
+      enddo
+
+c     write results, plot later
+      do i=0,ncells-1
+!        print *,(x(i)+x(i+1))*0.5d0,u(i)
+      enddo
+     
+      stop
+      end

examples/scalar_law/PROGRAM1/fortran/linearad.f

@@ -0,0 +1,68 @@
+c***********************************************************************
+c  Copyright 2006 John A. Trangenstein
+c
+c  This software is made available for research and instructional use 
+c  only. 
+c  You may copy and use this software without charge for these 
+c  non-commercial purposes, provided that the copyright notice and 
+c  associated text is reproduced on all copies.  
+c  For all other uses (including distribution of modified versions), 
+c  please contact the author at
+c    John A. Trangenstein
+c    Department of Mathematics
+c    Duke University
+c    Durham, NC 27708-0320
+c    USA
+c  or
+c    [email protected]
+c  
+c  This software is made available "as is" without any assurance that it
+c  is completely correct, or that it will work for your purposes.  
+c  Use the software at your own risk.
+c***********************************************************************
+      subroutine fluxderv(fc,lc,ifirst,ilast,
+     &  conservd,
+     &  dfdu    )
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      include "const.i"
+      include "linearad.i"
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      integer fc,lc,ifirst,ilast
+      double precision conservd(fc:lc)
+c
+      double precision dfdu(fc:lc)
+c
+      integer ic
+c     ******************************************************************
+c     compute derivative of flux with respect to conserved variable
+c     ******************************************************************
+      do ic=ifirst,ilast
+        dfdu(ic)=velocity
+      enddo
+
+      return
+      end 
+
+      subroutine riemann(fs,ls,ifirst,ilast,
+     &  left,right,
+     &  flux)
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      include "linearad.i"
+c+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+      integer fs,ls,ifirst,ilast
+c
+      double precision  
+     &  left(fs:ls+1),
+     &  right(fs:ls+1)
+c
+      double precision flux(fs:ls+1)
+      integer ie
+c     ******************************************************************
+c     the following riemann solver is hard-wired for linear advection
+c     with velocity > 0
+c     ******************************************************************
+      do ie=ifirst,ilast+1
+        flux(ie)=velocity*left(ie)
+      enddo
+      return
+      end 

examples/scalar_law/PROGRAM1/fortran/linearad.i

@@ -0,0 +1,26 @@
+c***********************************************************************
+c  Copyright 2006 John A. Trangenstein
+c
+c  This software is made available for research and instructional use 
+c  only. 
+c  You may copy and use this software without charge for these 
+c  non-commercial purposes, provided that the copyright notice and 
+c  associated text is reproduced on all copies.  
+c  For all other uses (including distribution of modified versions), 
+c  please contact the author at
+c    John A. Trangenstein
+c    Department of Mathematics
+c    Duke University
+c    Durham, NC 27708-0320
+c    USA
+c  or
+c    [email protected]
+c  
+c  This software is made available "as is" without any assurance that it
+c  is completely correct, or that it will work for your purposes.  
+c  Use the software at your own risk.
+c***********************************************************************
+      double precision jump,x_left,x_right,statelft,statergt,velocity,
+     &  dt_in
+      common/linearad/jump,x_left,x_right,statelft,statergt,velocity,
+     &  dt_in