Scrambler  1
IsoMHDWaves/problem.f90
Go to the documentation of this file.
00001 !#########################################################################
00002 !               
00003 !    Copyright (C) 2003-2012 Department of Physics and Astronomy,
00004 !                            University of Rochester,
00005 !                            Rochester, NY
00006 !
00007 !    problem.f90 of module IsoMHDWaves is part of AstroBEAR.
00008 !
00009 !    AstroBEAR is free software: you can redistribute it and/or modify    
00010 !    it under the terms of the GNU General Public License as published by 
00011 !    the Free Software Foundation, either version 3 of the License, or    
00012 !    (at your option) any later version.
00013 !
00014 !    AstroBEAR is distributed in the hope that it will be useful, 
00015 !    but WITHOUT ANY WARRANTY; without even the implied warranty of
00016 !    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00017 !    GNU General Public License for more details.
00018 !
00019 !    You should have received a copy of the GNU General Public License
00020 !    along with AstroBEAR.  If not, see <http://www.gnu.org/licenses/>.
00021 !
00022 !#########################################################################
00025 
00028 
00032 
00035 MODULE Problem
00036   USE DataDeclarations
00037   USE SplitRegions
00038   USE Shapes
00039   USE EOS
00040   USE Ambients
00041   USE RiemannSolvers
00042   USE Totals
00043   USE Fields
00044   IMPLICIT NONE
00045   SAVE
00046 
00047   PUBLIC ProblemModuleInit, ProblemGridInit, &
00048        ProblemBeforeStep, ProblemAfterStep, ProblemSetErrFlag, ProblemBeforeGlobalStep
00049   PRIVATE
00050   INTEGER, DIMENSIOn(:), ALLOCATABLE :: map
00051   REAL(KIND=qPREC), DIMENSION(:), ALLOCATABLE :: dq, q, w, dw
00052   REAL(KIND=qPREC), DIMENSION(:,:), ALLOCATABLE :: lambda, A, Lamb, Am, An, qExact
00053   REAL(KIND=qPREC), DIMENSION(:,:,:), ALLOCATABLE :: l, r
00054   REAL(KIND=qPREC) :: wave_speed
00055   REAL(KIND=qPREC) :: rho, p, vx, vy, vz, Bx, By, Bz, WaveStrength, wavespeed
00056   INTEGER :: WaveFamily
00057   
00058 CONTAINS
00059 
00061    SUBROUTINE ProblemModuleInit()      
00062       NAMELIST /ProblemData/ rho, p, vx, vy, vz, Bx, By, Bz, WaveFamily, WaveStrength
00063       OPEN(UNIT=PROBLEM_DATA_HANDLE, FILE='problem.data', STATUS="OLD")
00064       READ(PROBLEM_DATA_HANDLE,NML=ProblemData)
00065 
00066    END SUBROUTINE ProblemModuleInit
00067 
00070    SUBROUTINE ProblemGridInit(Info)
00071       TYPE(InfoDef) :: Info
00072       INTEGER :: i,j
00073       LOGICAL :: lFirstTime=.true.
00074       LOGICAL :: request_eigens
00075       INTEGER :: n(3)
00076       REAL(KIND=qPREC) :: left_fact, right_fact
00077       SAVE
00078       IF (lFirstTime) THEN
00079          lFirstTime=.false.
00080          IF (MPI_ID == 0) write(*,*) NrWaves, NrHydroVars
00081          ALLOCATE(w(1:NrHydroVars), dw(1:NrHydroVars), q(1:NrHydroVars),dq(1:NrHydroVars))
00082          ALLOCATE (r(nDim,NrWaves,NrWaves),l(nDim,NrWaves,NrWaves),lambda(nDim,NrWaves))
00083          ALLOCATE (A(NrHydroVars, NrHydroVars), Lamb(NrWaves, NrWaves), Am(NrWaves, NrWaves),An(NrWaves, NrWaves))
00084          w(1)=rho
00085          w(ivx)=vx
00086          IF (ivy /= 0) w(ivy)=vy
00087          IF (ivz /= 0) w(ivz)=vz
00088          IF (iE /= 0) w(iE)=P
00089          IF (iBx /= 0) w(iBx)=Bx
00090          IF (iBy /= 0) w(iBy)=By
00091          IF (iBz /= 0) w(iBz)=Bz
00092          write(*,*) 'sounds speed = ', sqrt(gamma*P/rho)
00093          write(*,*) w
00094          CALL calc_eigens(request_eigens, w, (/.true.,.true.,.true./),lambda, n, l, r,1,1,1,0)
00095 
00096          ALLOCATE(map(NrWaves))
00097          DO i=1,NrWaves
00098             map(i)=onedx_i(i)
00099             IF (lMHD .and. onedx_i(i) > iBx) map(i)=map(i)-1
00100          END DO
00101          write(*,*) 'map = ', map
00102          wavespeed=lambda(1,WaveFamily)
00103          write(*,*) 'wavespeeds', lambda(1,:)
00104          write(*,*) 'wave speed=', wavespeed
00105          write(*,*) 'loop should take ', (GxBounds(1,2)-GxBounds(1,1))/wavespeed
00106 !         left_fact=rho/sqrt(2d0*gamma*P/rho) / 0.1118033988749895
00107 !         right_fact=1d0/left_fact
00108 !         l=l*left_fact
00109 !         r=r*right_fact
00110          Lamb=0
00111          DO i=1,n(1)
00112             lamb(i,i)=lambda(1,i)
00113          END DO
00114 
00115          IF (MPI_ID == 0) THEN
00116             write(*,*) 'found ', n(1), 'waves'
00117            
00118             write(*,*) 'left eigen vectors are rows of'
00119             CALL OutputDoubleArray(transpose(l(1,1:n(1),map(:))))
00120             
00121             
00122             write(*,*) 'right eigen vectors are columns of'
00123             CALL OutputDoubleArray(r(1,1:n(1),map(:)))
00124             
00125             write(*,*) 'Check that L x R = I'
00126             CALL OutputDoubleArray(matmul(transpose(l(1,1:n(1),:)), r(1,1:n(1),:)))
00127             
00128 
00129 
00130             ! L x A x R = Lam
00131             ! L(-1) x Lam x R(-1) = R x Lam x L
00132             Am=matmul(matmul(transpose(l(1,1:n(1),:)), lamb(1:n(1),1:n(1))),r(1,1:n(1),:))
00133             write(*,*) 'Reconstruct matrix A'
00134             CALL OutputDoubleArray(Am)
00135             An=Am
00136             DO i=1,size(OneDx_i)
00137                An(:,map(i))=Am(:,i)
00138             END DO
00139             DO i=1,size(OneDx_i)
00140                Am(map(i),:)=An(i,:)
00141             END DO
00142             write(*,*) 'Reconstruct matrix A in normal index'
00143             CALL OutputDoubleArray(Am)
00144          END IF
00145 !         STOP
00146       END IF
00147       dw=0d0
00148       DO j=1,NrWaves
00149          dw(oneDx_i(j))=r(1,WaveFamily, j)
00150       END DO
00151       write(*,*) 'dw=', dw
00152 !      STOP
00153       dq=w+dw*1d-6
00154       CALL prim_to_cons(dq)
00155       q=w
00156       CALL prim_to_cons(q)
00157       write(*,*) 'dq=', (dq-q)/1d-6
00158 !      STOP
00159       DO i=1, Info%mX(1)
00160          Info%q(i,1,1,1:size(w)) = w+dw*WaveStrength*sin(2d0*Pi*real(i)/real(Info%mX(1)))
00161          CALL prim_to_cons(Info%q(i,1,1,:))
00162       END DO
00163       IF (MaintainAuxArrays) THEN
00164          IF (Info%mx(1) > Info%mX(2)) Info%aux(1:Info%mX(1),1:Info%mx(2)+1,1,2) = SPREAD(Info%q(1:Info%mX(1),1,1,iBy),2,Info%mX(2)+1)
00165          IF (Info%mX(2) > Info%mX(1)) Info%aux(1:Info%mX(1)+1,1:Info%mx(2),1,1) = SPREAD(Info%q(1,1:Info%mX(2),1,iBx),1,Info%mX(1)+1)
00166       END IF
00167 !      STOP
00168    END SUBROUTINE ProblemGridInit
00169 
00172    SUBROUTINE ProblemBeforeStep(Info)
00173       TYPE(InfoDef) :: Info
00174       INTEGER :: rmbc,i,j
00175       REAL(KIND=qPREC), DIMENSION(:,:,:), POINTER :: tempaux
00176       rmbc=levels(Info%level)%gmbc(1)
00177       IF (lMHD .AND. nDim == 2 .AND. Info%mX(1) == Info%mX(2)) THEN !Angled aux fields -  need to calculate potential
00178          ALLOCATE(tempaux(1:Info%mX(1)+1,1:Info%mX(2)+1,2))
00179          tempaux=Info%aux(1:Info%mX(1)+1,1:Info%mX(2)+1,1,1:2)
00180 
00181          DO j=-rmbc, Info%mX(2)+rmbc
00182             Info%aux(1-rmbc:j,Info%mx(2)+1-j,1,1)=Info%q(1,1,1,iBx)
00183             Info%aux(j+1:Info%mX(1)+1+rmbc,Info%mx(2)+1-j,1,1)=Info%q(Info%mX(1),Info%mX(2),1,iBx)
00184          END DO
00185          DO i=1-rmbc, Info%mX(1)+rmbc
00186             Info%aux(info%mx(1)+1-i,1-rmbc:i,1,2)=Info%q(1,1,1,iBy)
00187             Info%aux(info%mx(1)+1-i,i+1:Info%mX(2)+1+rmbc,1,2)=Info%q(Info%mX(1),Info%mX(2),1,iBy)
00188          END DO
00189          Info%aux(1:Info%mX(1)+1,1:Info%mX(2),1,1)=tempaux(1:Info%mX(1)+1,1:Info%mX(2),1)
00190          Info%aux(1:Info%mX(1),1:Info%mX(2)+1,1,2)=tempaux(1:Info%mX(1),1:Info%mX(2)+1,2)
00191 
00192          Info%q(1-rmbc:Info%mX(1)+rmbc,1-rmbc:Info%mX(2)+rmbc,1,iBx)=half*(Info%aux(1-rmbc:Info%mX(1)+rmbc,1-rmbc:Info%mX(2)+rmbc,1,1)+Info%aux(2-rmbc:Info%mX(1)+1+rmbc,1-rmbc:Info%mX(2)+rmbc,1,1))
00193          Info%q(1-rmbc:Info%mX(1)+rmbc,1-rmbc:Info%mX(2)+rmbc,1,iBy)=half*(Info%aux(1-rmbc:Info%mX(1)+rmbc,1-rmbc:Info%mX(2)+rmbc,1,2)+Info%aux(1-rmbc:Info%mX(1)+rmbc,2-rmbc:Info%mX(2)+1+rmbc,1,2))
00194          DEALLOCATE(tempaux)
00195       END IF
00196 
00197    END SUBROUTINE ProblemBeforeStep
00198 
00201    SUBROUTINE ProblemAfterStep(Info)
00202       TYPE(InfoDef) :: Info
00203       INTEGER :: i
00204       REAL(KIND=qPREC), DIMENSION(:), ALLOCATABLE :: err
00205       REAL(KIND=qPREC) :: phase
00206       ALLOCATE(err(1:NrHydroVars))
00207       IF (levels(info%level)%tnow+levels(info%level)%dt == final_time) THEN
00208          phase=-(2d0*Pi*(final_time*wavespeed)/(GxBounds(1,2)-GxBounds(1,1)))
00209          write(*,*) 'phase = ', phase
00210          OPEN(UNIT=11, FILE='out/data.curve', status='unknown')
00211          write(11,*) '# rho'         
00212          DO i=1, Info%mX(1)
00213             IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Info%q(i,1,1,1)
00214             IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Info%q(i,i,1,1)
00215          END DO
00216          write(11,*) 
00217          write(11,*) 
00218          write(11,*) '# rho_Exact'         
00219          DO i=1, Info%mX(1)
00220             IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(1)+dw(1)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00221             IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(1)+dw(1)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00222          END DO
00223          write(11,*) 
00224          write(11,*) 
00225 
00226 
00227          write(11,*) '# vx'
00228          DO i=1, Info%mX(1)
00229             IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Info%q(i,1,1,ivx)/Info%q(i,1,1,1)
00230             IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Info%q(i,i,1,ivx)/Info%q(i,i,1,1)
00231          END DO
00232          write(11,*) 
00233          write(11,*) 
00234 
00235          write(11,*) '# vx_Exact'         
00236          DO i=1, Info%mX(1)
00237             IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(ivx)+dw(ivx)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00238             IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(ivx)+dw(ivx)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00239          END DO
00240          write(11,*) 
00241          write(11,*) 
00242 
00243          write(11,*) '# P'
00244          DO i=1, Info%mX(1)
00245             IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, (Press(Info%q(i,1,1,:)))
00246             IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), (Press(Info%q(i,i,1,:)))
00247          END DO
00248          write(11,*) 
00249          write(11,*) 
00250          write(11,*) '# P_Exact'         
00251          IF (iE /= 0) THEN
00252             DO i=1, Info%mX(1)
00253                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(iE)+dw(iE)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00254                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(iE)+dw(iE)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00255             END DO
00256          ELSE
00257             DO i=1, Info%mX(1)
00258                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Iso_Speed2*(w(1)+dw(1)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1))))
00259                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Iso_Speed2*(w(1)+dw(1)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1))))
00260             END DO
00261          END IF
00262          write(11,*) 
00263          write(11,*) 
00264 
00265 
00266          IF (ivy /= 0) THEN
00267             write(11,*) '# vy'
00268             DO i=1, Info%mX(1)
00269                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Info%q(i,1,1,ivy)/Info%q(i,1,1,1)
00270                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Info%q(i,i,1,ivy)/Info%q(i,i,1,1)
00271             END DO
00272             write(11,*) 
00273             write(11,*) 
00274             write(11,*) '# vy_Exact'         
00275             DO i=1, Info%mX(1)
00276                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(ivy)+dw(ivy)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00277                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(ivy)+dw(ivy)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00278             END DO
00279             write(11,*) 
00280             write(11,*) 
00281          END IF
00282          IF (ivz /= 0) THEN
00283             write(11,*) '# vz'
00284             DO i=1, Info%mX(1)
00285                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Info%q(i,1,1,ivz)/Info%q(i,1,1,1)
00286                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Info%q(i,i,1,ivz)/Info%q(i,i,1,1)
00287             END DO
00288             write(11,*) 
00289             write(11,*) 
00290             write(11,*) '# vz_Exact'         
00291             DO i=1, Info%mX(1)
00292                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(ivz)+dw(ivz)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00293                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(ivz)+dw(ivz)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00294             END DO
00295             write(11,*) 
00296             write(11,*) 
00297 
00298          END IF
00299          IF (iBx /= 0) THEN
00300             write(11,*) '# Bx'
00301             DO i=1, Info%mX(1)
00302                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Info%q(i,1,1,iBx)
00303                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Info%q(i,i,1,iBx)
00304             END DO
00305             write(11,*) 
00306             write(11,*) 
00307             write(11,*) '# Bx_Exact'         
00308             DO i=1, Info%mX(1)
00309                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(iBx)+dw(iBx)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00310                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(iBx)+dw(iBx)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00311             END DO
00312             write(11,*) 
00313             write(11,*) 
00314 
00315          END IF
00316          IF (iBy /= 0) THEN
00317             write(11,*) '# By'
00318             DO i=1, Info%mX(1)
00319                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Info%q(i,1,1,iBy)
00320                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Info%q(i,i,1,iBy)
00321             END DO
00322             write(11,*) 
00323             write(11,*) 
00324             write(11,*) '# iBy_Exact'         
00325             DO i=1, Info%mX(1)
00326                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(iBy)+dw(iBy)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00327                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(iBy)+dw(iBy)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00328             END DO
00329             write(11,*) 
00330             write(11,*) 
00331 
00332          END IF
00333          IF (iBz /= 0) THEN
00334             write(11,*) '# Bz'
00335             DO i=1, Info%mX(1)
00336                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, Info%q(i,1,1,iBz)
00337                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), Info%q(i,i,1,iBz)
00338             END DO
00339             write(11,*) 
00340             write(11,*) 
00341             write(11,*) '# iBz_Exact'         
00342             DO i=1, Info%mX(1)
00343                IF (Info%mX(2) /= Info%mX(1)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx, w(iBz)+dw(iBz)*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1)))
00344                IF (Info%mX(1)==Info%mX(2)) write(11,*) Info%xbounds(1,1)+(REAL(i)-half)*levels(Info%level)%dx*sqrt(2d0), w(iBz)+dw(iBz)*WaveStrength*sin(phase+2d0*Pi*sqrt(2d0)*real(i)/real(Info%mX(1)))
00345             END DO
00346             write(11,*) 
00347             write(11,*) 
00348 
00349          END IF
00350          err=0d0
00351          DO i=1,Info%mX(1)
00352             CALL cons_to_prim(Info%q(i,1,1,:))
00353             err=err+abs(w+dw*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1))) - Info%q(i,1,1,1:NrHydroVars))
00354 !            write(*,*) i,abs(w+dw*WaveStrength*sin(phase+2d0*Pi*real(i)/real(Info%mX(1))) - Info%q(i,1,1,1:NrHydroVars))
00355          END DO
00356          err=err/real(Info%mX(1))
00357          write(*,*) 'dx, error=', levels(Info%level)%dx, sum(err)
00358       END IF
00359    END SUBROUTINE ProblemAfterStep
00360 
00363    SUBROUTINE ProblemSetErrFlag(Info)
00364       TYPE(InfoDef) :: Info
00365    END SUBROUTINE ProblemSetErrFlag
00366 
00367    SUBROUTINE ProblemBeforeGlobalStep(n)
00368       INTEGER :: n
00369    END SUBROUTINE ProblemBeforeGlobalStep
00370 
00371 END MODULE Problem
00372 
 All Classes Files Functions Variables