problem with smokefoil example

[Brian]

Hello,

I am having troubles running the examples:

amrcp latex/Chp2/method.3

and

amrcp latex/Chp2/method.4

I want to incorporate the smokefoil tracing into my own program but
when I run these example scripts I get the following error:

Error at line 155 of file /home/bmax/amrita/test/run_smokefoil:
unknown function 'W::VORT'!

Line 155 is:
minmax W::VORT[] -> min, max

error near:
[] -> min, max


It seems there is a problem when the script attempts to get the
vorticity when it reads the SootTraceImage procedure (for making
postscript files).

The version of Amrita I am running is v3.00 and here is the entire
example script I am running (latex/Chp2/method.3):


Latex2eHead {
dir = method
file = smokefoil.tex
}
fold::latex { title
fold::latex { top matter
\pagestyle{empty}
}
\centerline{\bf\Huge Shock Reflection From A Cylinder}
\bigskip
}
fold::amrita { run simulation
fold::amrita { prelims
EulerEquations {
add2W = VORT
}
fold::amrita { init resources
fold::amrita { get roe_fl
if(!&amrso("code/roe_fl")) then
fold::amrita { build code
fold::amrita { generate auxcode
fold::print { soot_trace
fold>file=code/soot_trace,recursion off,dollar
off
fold::amrita { daisy chain Ls
AddIntegrationStep {
step = soot (aka
Ls) :record_vort(IM,JM,NG,DX,DY,W)
}
set amrita:bcg::integration .= *Ls
}
fold::print { record_vort
*
*
* RECORD_VORT:
*
*
SUBROUTINE
RECORD_VORT(IM,JM,NG,DX,DY,W)
#include "$amrita:bcg::commonblk"
AMRINT IM,JM,NG
AMRDBL DX,DY
AMRDBL W(NW,amrCpatch(IM,JM,NG))
DO I=1,IM
DO J=1,JM
DUDY = (W(2,I,J+1)/W(1,I,J+1)
X -W(2,I,J-1)/W(1,I,J-1))/
(2*DY)
DVDX = (W(3,I+1,J)/W(1,I+1,J)
X -W(3,I-1,J)/W(1,I-1,J))/
(2*DX)
VORT = DVDX-DUDY
IF(ABS(VORT).GE.ABS(W(5,I,J)))
THEN
W(5,I,J) = ABS(VORT)
ENDIF
END DO
END DO
RETURN
END
}
}
}
BasicCodeGenerator {
solver = roe_fl
scheme = flux-limited'operator-split
auxcode= code/soot_trace
}
}
endif
}
fold::amrita { define DrawCylinderBody
proc DrawCylinderBody {
fold::amrita { parameters
xo = 60
yo = 0
r = 25
arc = 0 to 180
}
}
AmritaBlue
filled circle $xo,$yo,$r sector $arc
end proc
}
fold::amrita { define CylinderProblem
proc CylinderProblem {
fold::amrita { parameters
Ms = 1.25 # shock strength
Xs = 10.00 # shock position
lmax = 2 # grid levels
r = 2 # refinement ratio
rad = 25 # radius of cylinder
xo = 60 # centre of
yo = 0 # cylinder
}
}
fold::amrita { def Domain
def Domain
lscale 1
patch <+,1,w60,h60>
patch <+,1,w60,h60>
end def
}
fold::amrita { def BoundaryConditions
W'quiescent ::= <RHO=1,U=0,V=0,P=1>
fold::amrita { get post-shock state
ShockWave Ms=$Ms, state1=quiescent,\
state2=post_shock
}
def BoundaryConditions
Nbdy domain: reflect
Ebdy domain: extrapolate
Wbdy domain: prescribe W'post_shock
Sbdy domain: reflect
AXS::InternalBoundary {
fold::amrita { define internal boundary
PSI = ($rad-sqrt((X[]-($xo))**2+(Y[]-
($yo))**2))
}
}
end def
}
fold::amrita { def SolutionField
def SolutionField
setfield W'quiescent
setfield W'post_shock X[]<$Xs
end def
}
fold::amrita { def Mesh Adaption
def MeshAdaption
adaption on
lmax $lmax
r $r
end def
}
fold::amrita { def RefinementCriteria
def RefinementCriteria
DensityGradient
setflags [ooo|oxo|ooo] ON_BODY[]
end def
}
fold::amrita { now run it!
makefield
do l=1,$lmax
adapt
makefield
end do
}
end proc
}
fold::amrita { define SootTraceImage
proc SootTraceImage {
fold::amrita { parameters
SOOT = W::VORT # what to plot
exposure [0:1] = 1 # darkness of image
grid = {G} # select grid
tol = 0.2 # threshold for amp1
amp2
amp1 = 0.5 # weak features
amp2 = 0.25 # strong features
}
}
minmax $SOOT[] -> min, max
wt ::= abs(($SOOT[]-$min)/($max-$min))
greyshading ::= wt[]<$tol? wt[]**$amp1 : wt[]**$amp2
plot image $grid m<$exposure*greyshading[]>
end proc
}
}
plugin amr_sol
CylinderProblem {
fold::amrita { parameters
Ms = 2
Xs = 10
lmax = 1
r = 2
}
}
fold::amrita { init plotting
fold::amrita { force refinement
def RefinementCriteria
setflags [ooo|oxo|ooo] 1
end def
}
logfile logs/roe_fl
solver code/roe_fl
autoscale
postscript on
pushcwd $amrita:latex::PS
}
}
do phase=0,4
march 30 steps with cfl=0.8
fold::amrita { schlieren
plotfile sch$phase.ps
PlotDomain
SchlierenImage
DrawCylinderBody
plotfile
}
fold::amrita { soot trace
plotfile soot$phase.ps
PlotDomain
SootTraceImage
DrawCylinderBody
plotfile
}
end do
}
fold::latex { compose figures
fold::latex { schlierens
\noindent\huge schlieren images
\begin{center}
\begin{tabular}{ll}
fold::amrita { get fig
set fig = sch
}
fold::latex { snapshot t1
\begin{tabular}{l}
&latex::fig([width=7cm]{*fig`1})\\[5pt]
\LARGE snapshot t1
\end{tabular}
}
&
fold::latex { snapshot t2
\begin{tabular}{l}
&latex::fig([width=7cm]{*fig`2})\\[5pt]
\LARGE snapshot t2
\end{tabular}
}
\\[10pt]
fold::latex { snapshot t3
\begin{tabular}{l}
&latex::fig([width=7cm]{*fig`3})\\[5pt]
\LARGE snapshot t3
\end{tabular}
}
&
fold::latex { snapshot t4
\begin{tabular}{l}
&latex::fig([width=7cm]{*fig`4})\\[5pt]
\LARGE snapshot t4
\end{tabular}
}
\end{tabular}
\end{center}
}
fold::latex { soot trace
fold::amrita { get fig
set fig = soot4
}
\huge smoke foil
records vorticity history t1 through t4 \hfill
\begin{center}
\begin{tabular}{l}
&latex::fig([width=14cm]{*fig}) \\[5pt]
\end{tabular}
\end{center}
}
}
LatexTail
Latex

[Brian]

Also I am wondering if anyone has written any scripts, or has any
examples on plotting RHO*DZ/DT for the reactive euler equations case.

Thanks,

Brian

[James Quirk]

Brian,

On Mon, 16 Jun 2008, Brian wrote:

>
> Hello,
>
> I am having troubles running the examples:
>
> amrcp latex/Chp2/method.3
>
> and
>
> amrcp latex/Chp2/method.4
>
> I want to incorporate the smokefoil tracing into my own program but
> when I run these example scripts I get the following error:
>
> Error at line 155 of file /home/bmax/amrita/test/run_smokefoil:
> unknown function 'W::VORT'!
>
> Line 155 is:
> minmax W::VORT[] -> min, max
>
> error near:
> [] -> min, max
>

The problem is caused by a bug in the machinery that sets up an
EquationSet. For the fix, see my posting to Matei in the thread
"shock/detonation interaction with a ramp", May 27.

James

[James Quirk]

Brian,

On Mon, 16 Jun 2008, Brian wrote:

>
> Also I am wondering if anyone has written any scripts, or has any
> examples on plotting RHO*DZ/DT for the reactive euler equations case.

You'd have to elaborate a bit more on exactly what you want to plot.
Because of the way the mesh adaption works, it's not possible to provide
a general mechanism for accessing time derivatives in the flow solution.
But if you tell me what you want to do, it should be possible to introduce
some custom code to do the trick.

James
>
> Thanks,
>
> Brian
>

[Brian]

Thanks James,

I fixed the bug and the script works fine now. Do you have any
examples of using the Add2W = VORT for the reactive euler case? I
have tried but keeping *** MALLOC FAILED IN GRAB_DSTR *** error... :(

The problem I am working on is examining the stability of a detonation
wave after it passes over a cylinder (or a couple cylinders). It
would be nice to trace the regions of maximum vorticity as in the
example scripts as well as regions of high reaction rates (RHO*DZ/DT)
integrated over time.

Brian

[James Quirk]

Brian,

> Thanks James,
>
> I fixed the bug and the script works fine now. Do you have any
> examples of using the Add2W = VORT for the reactive euler case? I

No, but Gary does. But you will see its use is exactly the
same as the for the plain EulerEquations.

> have tried but keeping *** MALLOC FAILED IN GRAB_DSTR *** error... :(

It sounds like you will need to invoke the routine ArraySizes
to increase Amr_sol's storage allocation. For example, you could try:

ArraySizes {
NGIxJ = 1000000
NGD = 4000
}

>
> The problem I am working on is examining the stability of a detonation
> wave after it passes over a cylinder (or a couple cylinders). It
> would be nice to trace the regions of maximum vorticity as in the
> example scripts as well as regions of high reaction rates (RHO*DZ/DT)
> integrated over time.

You could keep track of the rate term in exactly the same way that
the vorticity is tracked. But rather than me spoon-feed you,
I would prefer for you to have a stab on your own, first.

James

[Brian]

Thanks James,

I will try that. I also have a problem where I may want to convect
some inert scalar, i'll call it alpha, throughout the domain so that I
can track one gas expanding into another. For example I will
initially set gas A to have a value of alpha=0 and gas B will have a
value of alpha=1. Is there a way I can add this scalar to the
equation set for the EulerEquations (not reactive) and take advantage
of the built in solvers (Roe / Godonov), or do I have to write my own
discretization for convecting alpha?

Thanks for your help!

Brian

[James Quirk]

Brian,

On Sun, 22 Jun 2008, Brian wrote:

>
> Thanks James,
>
> I will try that. I also have a problem where I may want to convect
> some inert scalar, i'll call it alpha, throughout the domain so that I
> can track one gas expanding into another. For example I will
> initially set gas A to have a value of alpha=0 and gas B will have a
> value of alpha=1. Is there a way I can add this scalar to the
> equation set for the EulerEquations (not reactive) and take advantage
> of the built in solvers (Roe / Godonov), or do I have to write my own
> discretization for convecting alpha?

BCG could be coerced into doing what you want, but as I've
just pointed out to Gary: BCG's internal construction is
not in the best of shapes and so it would probably be more
trouble than it's worth. Here I should point out that
BCG really is "basic" and it was never intended to be all
things to all people.

Now I understand from Matei that you are signed up to do a Ph.D.
and so if you intend to be doing reactive flow simulations
for any length of time, you should really take the time
to dissect some of the solver's that BCG constructs, so
as to craft your own custom solver. The amount of code
is actually quite small and anything you produce
can always be absorbed back into BCG. And it would be
nice to think that you wanted to go the extra mile,
with the aim of leaving a legacy for Matei's next student.

James

[Gary Sharpe]

Brian, also I recall the Charles was doing a similar thing - he needed
two advecting
scalars to track where a particle started. Have a look through the
newsgroup postings
- there may be some hints.