Hello,
I am doing a stability calculation for different profiles.
The test case that I have used is Gaussian, where the velocity
and gradient of PV is defined using something like this,
Ub = Function(V).interpolate(Uj*exp(-(x[0]-Ly/2)**2))
dQb =
Function(V).interpolate(-2.*Uj*exp(-(x[0]-Ly/2)**2)*(2.*(x[0]-Ly/2)**2
- 1.))
I have compared this with the results of using a spectral method
and the results agree very well.
I am now considering a less trivial example using something
called a Bump function. Yes, technical name is Bump. Thanks to
previous help I have figured out how to define it piecewise
(thanks again).
Ub = Function(V).interpolate(conditional(le(abs(x[0]-Ly/2.),
1.0), exp(-1./(1. -(x[0]-Ly/2)**2)), 0.0))
dQb =
Function(V).interpolate(conditional(le(abs(x[0]-Ly/2.), 1.0),
-2.*exp(-1./(1. -(x[0]-Ly/2)**2))*(3.*(x[0]-Ly/2.)**4 -
1)/((x[0]-Ly/2.)**2-1)**4, 0.0))
However, when I try this I get zero growth rate, which does
not agree with the spectral method.
Actually, in the spectral method we need to zero out the Nans
and Infs in expressions, and that works out nicely. Is there
an easy way to convert Nans and Infs to zero in interpolate?
I tried using both a Krylov method and Lapack. The first
gives zero, which I don't believe, and the second fails, I
presume because I am dividing by zero somewhere.
Any help would be greatly appreciated!
Cheers, Francis