Commit 77a96f60 authored by Stefan Kerkemeier's avatar Stefan Kerkemeier Committed by GitHub

Merge pull request #59 from stgeke/par_imprv

Fix for different ps solvers
parents 1d715892 bd391251
......@@ -197,13 +197,8 @@ void iniparser_dump(const dictionary * d, FILE * f)
if (d==NULL || f==NULL) return ;
for (i=0 ; i<d->size ; i++) {
if (d->key[i]==NULL)
continue ;
if (d->val[i]!=NULL) {
fprintf(f, "%s = [%s]\n", d->key[i], d->val[i]);
} else {
/*fprintf(f, "[%s]=UNDEF\n", d->key[i]);*/
}
if (d->key[i]==NULL) continue ;
fprintf(f, "%s = [%s]\n", d->key[i], d->val[i]);
}
return ;
}
......@@ -579,7 +574,8 @@ static line_status iniparser_line(
} else if (line[0]=='#' || line[0]==';') {
/* Comment line */
sta = LINE_COMMENT ;
} else if (line[0]=='[' && line[len-1]==']') {
/* } else if (line[0]=='[' && line[len-1]==']') { */
} else if (line[0]=='[') {
/* Section name */
sscanf(line, "[%[^]]", section);
strstrip(section);
......
integer cv_lysize
parameter (cv_lysize=lx1*ly1*lz1*lelt*ldimt+1)
integer cv_ystart, cv_1fld, cv_lfld
COMMON /ICVODE/ cv_ystart, cv_1fld, cv_lfld
integer cv_nfld, cv_nadd_usr
COMMON /ICVODE/ cv_nfld, cv_nadd_usr
integer*8 cv_nglobal
COMMON /ILCVODE/ cv_nglobal
......@@ -11,16 +11,16 @@
COMMON /LCVODE/ ifcvodeinit
real cv_atol,cv_rtol
real cv_dtlag,cv_abmsh,cv_ab,cv_dtnek
real cv_dtlag,cv_abmsh,cv_ab
real cv_time,cv_timel,cv_dtnek
COMMON /RCVODE/ cv_atol,cv_rtol,
& cv_dtlag(3),cv_abmsh(3),cv_ab(3),cv_dtnek
& cv_dtlag(3),cv_abmsh(3),cv_ab(3),
& cv_time,cv_timel,cv_dtnek
INTEGER IGSTYPE,METH,ITMETH,IATOL,ITASK
INTEGER IGSTYPE,ITMETH,IATOL,ITASK
PARAMETER(
& METH = 2, ! BDF
& ITMETH = 2, ! Newton iter
& IATOL = 1, ! 1: CV_SS / 2: CV_SV)
& IGSTYPE = 1, ! GS 1: modified 2: classical
& ITASK = 1 ! overshoot + interpolate in time
& IGSTYPE = 1 ! GS 1: modified 2: classical
& )
......@@ -22,7 +22,7 @@ C
$ ,NKTONV,NHIS,LOCHIS(4,lhis+maxobj)
$ ,IPSCAL,NPSCAL,IPSCO, ifldmhd
$ ,IRSTV,IRSTT,IRSTIM,NMEMBER(MAXOBJ),NOBJ
$ ,NGEOM
$ ,NGEOM,IDPSS(LDIMT)
C
COMMON /INPUT3/ IF3D,IFFLOW,IFHEAT,IFTRAN,IFAXIS,IFSTRS,IFSPLIT
$ ,IFMGRID,IFADVC(LDIMT1),IFTMSH(0:LDIMT1)
......@@ -34,6 +34,7 @@ C
$ ,IFCVODE,IFLOMACH,IFEXPLVIS,IFSCHCLOB,IFUSERVP
$ ,IFCYCLIC,IFMOAB,IFCOUP, IFVCOUP, IFUSERMV,IFREGUO
$ ,IFXYO_,ifaziv,IFNEKNEK,IFNEKNEKM,ifdg
$ ,IFCVFLD(LDIMT1),IFDP0DT
c
LOGICAL IF3D,IFFLOW,IFHEAT,IFTRAN,IFAXIS,IFSTRS,IFSPLIT
$ ,IFMGRID,IFADVC ,IFTMSH
......@@ -45,6 +46,7 @@ c
$ ,IFCVODE,IFLOMACH,IFEXPLVIS,IFSCHCLOB,IFUSERVP
$ ,IFCYCLIC,IFMOAB,IFCOUP, IFVCOUP, IFUSERMV,IFREGUO
$ ,IFXYO_,ifaziv,IFNEKNEK,IFNEKNEKM,ifdg
$ ,IFCVFLD,IFDP0DT
LOGICAL IFNAV
EQUIVALENCE (IFNAV, IFADVC(1))
......
......@@ -3,7 +3,7 @@ c Define all valid keys for .par file here
c Note: Keys have to be in captial letters
c
integer PARDICT_NKEYS
parameter(PARDICT_NKEYS = 68)
parameter(PARDICT_NKEYS = 75)
character*132 pardictkey(PARDICT_NKEYS)
data
......@@ -68,10 +68,17 @@ c
& pardictkey(59) / 'GENERAL:USERPARAM09' /
& pardictkey(60) / 'GENERAL:USERPARAM10' /
& pardictkey(61) / 'SCALAR%%' /
& pardictkey(62) / 'SCALAR%%:RHOCP' /
& pardictkey(63) / 'SCALAR%%:CONDUCTIVITY' /
& pardictkey(64) / 'SCALAR%%:WRITETOFIELDFILE' /
& pardictkey(65) / 'TEMPERATURE:SOLVER' /
& pardictkey(66) / 'TEMPERATURE:ABSOLUTETOL' /
& pardictkey(67) / 'TEMPERATURE:RELATIVETOL' /
& pardictkey(62) / 'SCALAR%%:SOLVER' /
& pardictkey(63) / 'SCALAR%%:RESIDUALTOL' /
& pardictkey(64) / 'SCALAR%%:RHO' /
& pardictkey(65) / 'SCALAR%%:DIFFUSIVITY' /
& pardictkey(66) / 'SCALAR%%:WRITETOFIELDFILE' /
& pardictkey(67) / 'TEMPERATURE:SOLVER' /
& pardictkey(68) / 'TEMPERATURE:RESIDUALTOL' /
& pardictkey(69) / 'PROBLEMTYPE:DP0DT' /
& pardictkey(70) / 'CVODE' /
& pardictkey(71) / 'CVODE:ABSOLUTETOL' /
& pardictkey(72) / 'CVODE:RELATIVETOL' /
& pardictkey(73) / 'CVODE:STIFF' /
& pardictkey(74) / 'CVODE:MODE' /
& pardictkey(75) / 'CVODE:PRECONDITIONER' /
......@@ -121,6 +121,7 @@ c Fill up user defined forcing function and collocate will the
c mass matrix on the Gauss-Lobatto mesh.
include 'SIZE'
include 'INPUT'
include 'MASS'
include 'SOLN'
include 'TSTEP'
......@@ -129,9 +130,11 @@ c mass matrix on the Gauss-Lobatto mesh.
time = time-dt ! Set time to t^n-1 for user function
call rzero ( bq(1,1,1,1,ifield-1) , n)
call setqvol ( bq(1,1,1,1,ifield-1) )
call col2 ( bq(1,1,1,1,ifield-1) ,bm1,n)
call rzero (bq(1,1,1,1,ifield-1),n)
call setqvol (bq(1,1,1,1,ifield-1))
if (ifcvfld(ifield)) call copy (tlag(1,1,1,1,1,ifield-1),
& bq (1,1,1,1,ifield-1),n)
call col2 (bq(1,1,1,1,ifield-1) ,bm1,n)
time = time+dt ! Restore time
......
......@@ -289,7 +289,12 @@ C
do i=1,NPSCL2
IFTMSH(i) = .false.
IFADVC(i) = .false.
enddo
enddo
do i=1,NPSCL1
IDPSS(i) = 0 ! use Helmholtz for passive scalars
enddo
IFFLOW = .false.
IFHEAT = .false.
IFTRAN = .false.
......
subroutine add_fcvfun_usr(ydot,jstart)
subroutine fcvfun_usr(ydot,jstart)
c
c plug-in specific contributions to rhs
c user specific contributions to rhs
c
include 'SIZE'
include 'TOTAL'
real ydot(1)
ntotv = nx1*ny1*nz1*nelv
if (ifvcor .and. iflomach) then
dd = gamma0 ! CVref/CPref ! Note CVref denotes the inverse CPref
dd = (dd - 1.)/dd
xfacr= dd * dp0thdt
do i = 1,ntotv
ydot(i) = ydot(i) + xfacr/vtrans(i,1,1,1,2)
enddo
else
dp0thdt = 0.0
endif
if (ifvcor .and. iflomach) ydot(jstart) = dp0thdt
return
end
c----------------------------------------------------------------------
subroutine cv_unpack_sol(y)
subroutine cvunpack_usr(y,jstart)
c
c copy the cvode solution (y) back to the internal nek array (t)
c
......@@ -38,21 +20,10 @@ c
real y(1)
nxyz = nx1*ny1*nz1
j = 1
do i=cv_1fld,cv_lfld
ntot = nxyz*nelfld(i)
call copy(t(1,1,1,1,i-1),y(j),ntot)
j = j + ntot
enddo
if (ifvcor .and. iflomach) p0th = y(j)
return
end
c----------------------------------------------------------------------
subroutine cv_pack_sol(y)
subroutine cvpack_usr(y,jstart)
c
c copy the internal nek array (t) to the cvode solution (y)
c
......@@ -62,16 +33,5 @@ c
real y(1)
nxyz = nx1*ny1*nz1
j = 1
do i=cv_1fld,cv_lfld
ntot = nxyz*nelfld(i)
call copy (y(j),t(1,1,1,1,i-1),ntot)
j = j + ntot
enddo
if (ifvcor .and. iflomach) y(j) = p0th
return
end
This diff is collapsed.
......@@ -91,12 +91,6 @@ c COMMON /SCRCG/ DUMM10(LX1,LY1,LZ1,LELT,1)
call bcmask ! Set BC masks for Dirichlet boundaries.
if (ifcvode.and.nsteps.gt.0) then
n_aux = 0
if(iflomach .and. ifvcor) n_aux = 1 ! Thermodynamic pressure
call cv_setsize(2,nfield,n_aux) ! Set size for CVODE solver
endif
if (fintim.ne.0.0.or.nsteps.ne.0)
$ call geneig(igeom) ! eigvals for tolerances
......@@ -118,20 +112,18 @@ c COMMON /SCRCG/ DUMM10(LX1,LY1,LZ1,LELT,1)
endif
endif
call init_plugin ! Initialize optional plugin
if(ifcvode) call cv_setsize(0)
if(nio.eq.0) write(6,*) 'call usrdat3'
call usrdat3
if(nio.eq.0) write(6,'(A,/)') ' done :: usrdat3'
call cmt_switch ! Check if compiled with cmt
if (ifcmt) then ! Initialize CMT branch
call nek_cmt_init
if (nio.eq.0) write(6,*)'Initialized DG machinery'
endif
call setics ! Set initial conditions
call setprop ! Compute field properties
......@@ -144,7 +136,7 @@ c COMMON /SCRCG/ DUMM10(LX1,LY1,LZ1,LELT,1)
if(nio.eq.0) write(6,'(A,/)') ' done :: userchk'
endif
if(ifcvode .and. nsteps.gt.0) call cv_init ! Initialize CVODE
if (ifcvode .and. nsteps.gt.0) call cv_init
call comment
call sstest (isss)
......@@ -233,6 +225,7 @@ c check for post-processing mode
RETURN
END
c-----------------------------------------------------------------------
subroutine nek_advance
......@@ -260,24 +253,24 @@ c-----------------------------------------------------------------------
! within cvode we use the lagged wx for
! extrapolation, that's why we have to call it before gengeom
if (ifheat .and. ifcvode) call heat (igeom)
if (ifheat .and. ifcvode) call heat_cvode (igeom)
if (ifgeom) then
call gengeom (igeom)
call geneig (igeom)
endif
if (ifheat .and. .not.ifcvode) call heat (igeom)
if (ifheat) call heat (igeom)
if (igeom.eq.2) then
call setprop
call qthermal(1)
call setprop
if (iflomach) call qthermal(.true.,.false.,dummy)
endif
if (ifflow) call fluid (igeom)
if (ifmvbd) call meshv (igeom)
if (param(103).gt.0) call q_filter (param(103))
call setup_convect (igeom) ! Save convective velocity _after_ filter
if (ifflow) call fluid (igeom)
if (ifmvbd) call meshv (igeom)
if (param(103).gt.0) call q_filter (param(103))
call setup_convect (igeom) ! Save convective velocity _after_ filter
enddo
else ! PN-2/PN-2 formulation
......@@ -318,6 +311,7 @@ c-----------------------------------------------------------------------
return
end
c-----------------------------------------------------------------------
subroutine nek_end
......
......@@ -82,6 +82,7 @@ C
NTOT=NX2*NY2*NZ2*LELT
CALL RZERO(USRDIV,NTOT)
CALL RZERO(QTL,NTOT)
RETURN
END
......@@ -240,9 +241,6 @@ C
call exitt
endif
if(abs(PARAM(16)).ge.2) IFCVODE = .true.
C
C Check accuracy requested.
C
IF (TOLREL.LE.0.) TOLREL = 0.01
......@@ -747,7 +745,7 @@ C-----------------------------------------------------------------------
ts = dnekclock()
if(nio.eq.0 .and. igeom.eq.2)
& write(6,'(13x,a)') 'Solving for fluid' !,ifsplit,iftran,ifnav
& write(*,'(13x,a)') 'Solving for fluid'
if (ifsplit) then
......@@ -784,8 +782,8 @@ c - Velocity/stress formulation.
endif
if(nio.eq.0 .and. igeom.ge.2)
& write(*,'(11x,a,11x,1p1e12.4,a)') ' Fluid done',
& dnekclock()-ts
& write(*,'(4x,i7,a,1p2e12.4)')
& istep,' Fluid done',time,dnekclock()-ts
return
end
......@@ -813,40 +811,49 @@ C
ts = dnekclock()
if (nio.eq.0 .and. igeom.eq.1)
& write(*,'(13x,a)') 'Solving for heat'
if (ifcvode) then
if(igeom.eq.2) call cdscal_cvode(igeom)
if (nio.eq.0 .and. igeom.eq.2)
& write(*,'(13x,a)') 'Solving for Hmholtz scalars'
elseif (ifsplit) then
do ifield=2,nfield
do ifield = 2,nfield
if (idpss(ifield-1).eq.0) then ! helmholtz
intype = -1
if (.not.iftmsh(ifield)) imesh = 1
if ( iftmsh(ifield)) imesh = 2
call unorm
call settolt
call cdscal (igeom)
enddo
call cdscal(igeom)
endif
enddo
else ! PN-PN-2
if (nio.eq.0 .and. igeom.eq.2)
& write(*,'(4x,i7,a,1p2e12.4)')
& istep,' Scalars done',time,dnekclock()-ts
do ifield=2,nfield
intype = -1
if (.not.iftmsh(ifield)) imesh = 1
if ( iftmsh(ifield)) imesh = 2
call unorm
call settolt
call cdscal (igeom)
enddo
return
end
c-----------------------------------------------------------------------
subroutine heat_cvode (igeom)
C
include 'SIZE'
include 'INPUT'
include 'TSTEP'
include 'TURBO'
include 'DEALIAS'
endif
real*8 ts, dnekclock
ts = dnekclock()
if (igeom.ne.2) return
if (nio.eq.0)
& write(*,'(13x,a)') 'Solving for CVODE scalars'
call cdscal_cvode(igeom)
if (nio.eq.0 .and. igeom.ge.2)
& write(*,'(11x,a,12x,1p1e12.4,a)') ' Heat done',
& dnekclock()-ts
if (nio.eq.0)
& write(*,'(4x,i7,a,1p2e12.4)')
& istep,' CVODE scalars done',time,dnekclock()-ts
return
end
......
subroutine init_plugin
return
end
......@@ -45,7 +45,7 @@ genbox.o gmres.o hsmg.o convect.o induct.o perturb.o \
navier5.o navier6.o navier7.o navier8.o fast3d.o fasts.o calcz.o \
byte.o chelpers.o byte_mpi.o postpro.o \
cvode_driver.o nek_comm.o \
init_plugin.o setprop.o qthermal.o cvode_aux.o makeq_aux.o \
vprops.o qthermal.o cvode_aux.o makeq_aux.o \
papi.o nek_in_situ.o \
readat_new.o finiparser.o iniparser.o dictionary.o
################################################################################
......@@ -226,8 +226,7 @@ $(OBJDIR)/gfdm_op.o :$S/gfdm_op.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/coef.o :$S/coef.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/plan4.o :$S/plan4.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/qthermal.o :$S/qthermal.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/setprop.o :$S/setprop.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/init_plugin.o :$S/init_plugin.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/vprops.o :$S/vprops.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/cvode_driver.o :$S/cvode_driver.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/cvode_aux.o :$S/cvode_aux.f; $(F77) -c $(FL2) $< -o $@
$(OBJDIR)/makeq.o :$S/makeq.f; $(F77) -c $(FL2) $< -o $@
......
......@@ -18,9 +18,6 @@ CC="mpicc"
# (set PPLIST=? to get a list of available symbols)
#PPLIST="?"
# plug-in list
#PLUGIN_LIST=""
# OPTIONAL SETTINGS
# -----------------
......
......@@ -378,7 +378,7 @@ if [ $? -eq 0 ]; then
fi
echo $PPLIST | grep 'BGQ' >/dev/null
if [ $? -eq 0 ]; then
MXM_USER="mxm_std.o mxm_bgq.o blas.o"
MXM_USER="mxm_std.o mxm_bgq.o "
OPT_FLAGS_STD="-O3"
OPT_FLAGS_MAG="-O3"
fi
......
......@@ -10,31 +10,33 @@ C !! NOTE: Do not change the content of the array BQ until the current
logical if_conv_std
common /SCRUZ/ w1(lx1,ly1,lz1,lelt)
if_conv_std = .true.
if (ifmhd.and.ifaxis) if_conv_std = .false. ! conv. treated in induct.f
nxyz = nx1*ny1*nz1
ntot = nxyz*nelv
if_conv_std = .true.
if (ifmhd.and.ifaxis) if_conv_std = .false. ! conv. treated in induct.f
! source terms
call makeq_aux ! nekuq, etc.
if(ifcvode .and. ifmvbd) then
if(ifcvfld(ifield) .and. ifmvbd) then
call sub2 (vx, wx, ntot)
call sub2 (vy, wy, ntot)
call sub2 (vz, wz, ntot)
endif
! update fine grid velocity
if (param(99).gt.0) call set_convect_new(vxd,vyd,vzd,vx,vy,vz)
! convection term
if (param(99).gt.0) call set_convect_new(vxd,vyd,vzd,vx,vy,vz)
if (ifadvc(ifield).and..not.ifchar.and.if_conv_std) call convab
if (ifcvode .and. ifmvbd) then
if (ifcvfld(ifield) .and. ifmvbd) then
call add2 (vx, wx, ntot) ! add back W to V for mvmesh
call add2 (vy, wy, ntot)
call add2 (vz, wz, ntot)
endif
if (iftran) then
if(ifcvode) then
if(ifcvfld(ifield)) then
ntot = nx1*ny1*nz1*nelfld(ifield)
call wlaplacian(w1,t(1,1,1,1,ifield-1),vdiff(1,1,1,1,ifield),
& ifield)
......
......@@ -8,11 +8,6 @@
if_conv_std = .true.
if (ifmhd.and.ifaxis) if_conv_std = .false. ! conv. treated in induct.f
if(ifcvode .and. ifield.eq.2) then
call setprop
call qthermal(0)
endif
call whatfld (ifturb)
if (ifturb) call maketq ! zero bq
if (.not.ifturb .and. if_conv_std) call makeuq !zero bq
......
......@@ -1887,4 +1887,28 @@ c
c
return
end
cc-----------------------------------------------------------------------
subroutine transpose(a,lda,b,ldb)
real a(lda,1),b(ldb,1)
c
do j=1,ldb
do i=1,lda
a(i,j) = b(j,i)
enddo
enddo
return
end
c-----------------------------------------------------------------------
subroutine transpose1(a,n)
real a(n,n)
c
do j=1,n
do i=j+1,n
ta = a(i,j)
a(i,j) = a(j,i)
a(j,i) = ta
enddo
enddo
return
end
c-----------------------------------------------------------------------
This diff is collapsed.
......@@ -29,28 +29,36 @@ c
#ifdef BGQ
if (n2 .eq. 8 .and. mod(n1,4) .eq. 0 ) then
if (n2 == 8 .and. mod(n1,4) == 0 &
c .and. MOD(LOC(a),tt)==0 &
c .and. MOD(LOC(b),tt)==0 &
c .and. MOD(LOC(c),tt)==0 &
) then
call mxm_bgq_8(a, n1, b, n2, c, n3)
return
endif
if (n2 .eq. 16 .and. mod(n1,4) .eq. 0 ) then
if (n2 == 16 .and. mod(n1,4) == 0 &
c .and. MOD(LOC(a),tt)==0 &
c .and. MOD(LOC(b),tt)==0 &
c .and. MOD(LOC(c),tt)==0 &
) then
call mxm_bgq_16(a, n1, b, n2, c, n3)
return
endif
tt = 32
if (n2 .eq. 10 .and. mod(n1,4) .eq. 0 .and. mod(n3,2) .eq. 0
+ .and. MOD(LOC(a),tt).eq.0
+ .and. MOD(LOC(b),tt).eq.0
+ .and. MOD(LOC(c),tt).eq.0
+ ) then
if (n2 == 10 .and. mod(n1,4) == 0 .and. mod(n3,2) == 0 &
.and. MOD(LOC(a),tt)==0 &
.and. MOD(LOC(b),tt)==0 &
.and. MOD(LOC(c),tt)==0 &
) then
call mxm_bgq_10(a, n1, b, n2, c, n3)
return
endif
if (n2 .eq. 6 .and. mod(n1,4) .eq. 0 .and. mod(n3,2) .eq. 0
+ .and. MOD(LOC(a),tt).eq.0
+ .and. MOD(LOC(b),tt).eq.0
+ .and. MOD(LOC(c),tt).eq.0
+ ) then
if (n2 == 6 .and. mod(n1,4) == 0 .and. mod(n3,2) == 0 &
.and. MOD(LOC(a),tt)==0 &
.and. MOD(LOC(b),tt)==0 &
.and. MOD(LOC(c),tt)==0 &
) then
call mxm_bgq_6(a, n1, b, n2, c, n3)
return
endif
......
......@@ -3414,17 +3414,6 @@ c
c
return
end
c-----------------------------------------------------------------------
subroutine transpose(a,lda,b,ldb)
real a(lda,1),b(ldb,1)
c
do j=1,ldb
do i=1,lda
a(i,j) = b(j,i)
enddo
enddo
return
end
c-----------------------------------------------------------------------
subroutine convop(conv,fi)
C
......
......@@ -117,7 +117,7 @@ c
$ call filterq(vzp(1,j),intv,nx1,nz1,wk1,wk2,intt,if3d,wmax)
ifield = 2
if (ifheat .and. .not.ifcvode)
if (ifheat .and. .not.ifcvfld(ifield))
$ call filterq(tp(1,j,1),intv,nx1,nz1,wk1,wk2,intt,if3d,wmax)
enddo
......@@ -1477,19 +1477,6 @@ c
c
return
end
c-----------------------------------------------------------------------
subroutine transpose1(a,n)
real a(n,n)
c
do j=1,n
do i=j+1,n
ta = a(i,j)
a(i,j) = a(j,i)
a(j,i) = ta
enddo
enddo
return
end
c-----------------------------------------------------------------------
subroutine get_exyz(ex,ey,ez,eg,nelx,nely,nelz)
integer ex,ey,ez,eg
......
......@@ -131,7 +131,7 @@ c Trigger history output only if prefix = 'his' pff 8/18/05
if (iiidmp.ne.0.or.lastep.eq.1.or.timdump.eq.1.) ifdoit=.true.
if (ifdoit.and.nio.eq.0)write(6,*)'call outfld: ifpsco:',ifpsco(1)
c if (ifdoit.and.nio.eq.0)write(6,*)'call outfld: ifpsco:',ifpsco(1)
if (ifdoit) call outfld(prefix)
call outhis(ifhis)
......
subroutine qthermal(iflag)
subroutine qthermal(ifvext,ifqvol,qvol)
C Compute the thermal divergence QTL
C
......@@ -14,6 +14,9 @@ C energy equation expressed in terms of temperature.
include 'SIZE'
include 'TOTAL'
logical ifvext,ifqvol
real qvol(1)
common /scrns/ w1(lx1,ly1,lz1,lelt)
$ ,w2(lx1,ly1,lz1,lelt)
$ ,w3(lx1,ly1,lz1,lelt)
......@@ -28,17 +31,17 @@ C energy equation expressed in terms of temperature.
nxyz = nx1*ny1*nz1
ntot = nxyz*nelv
if (.not.iflomach) then
call rzero(qtl,ntot)
return
endif