projecting with inhomogeneous boundary conditions about shenfun HOT 5 CLOSED

Interesting. A few things go wrong here and I don't think these things are very well documented.

1. P0 uses boundary conditions (-1, 1), but the function you are trying to approximate pr has boundary conditions (-sin(y), sin(y)). This can be easily fixed using

import sympy as sp
x, y = sp.symbols('x,y', real=True)
P0   = Basis(N[0], family='C', domain=(-L[0]/2,L[0]/2), bc=(a*sp.sin(y), b*sp.sin(y)), scaled=True)

2. The gradient grad(pk) has homogeneous boundary conditions for it's first item and inhomogeneous (a*sp.cos(y), b*sp.cos(y)) for its second item. I would suggest projecting the gradient to the regular Chebyshev or Legendre space that has no boundary conditions. Much more convenient:-)

P1   = Basis(N[0], family='C', domain=(-L[0]/2,L[0]/2))
TF1  = TensorProductSpace(comm, (P1, F1), axes=(0, 1))
TV1  = VectorTensorProductSpace(TF1) 
ur = project(grad(pk), TV1).backward()

The reason your code fails is, like you suggest, that one cannot project to a Dirichlet space with inhomogeneous boundary conditions. I think I could make that work quite easily, but I guess it has not come up until now because it has always been a better idea to project to the regular, non-restricted, spaces, since boundary conditions are not required for a projection.

from shenfun.

Ah, that makes so much sense. I am sorry for not figuring this out myself but what you say is very clear.

I don't think there is anything for you to change, it's more like the user, me, needing to be clearer as to what I want to do.

The P0 space is for my streamfunction and I was foolishly using this to build my space for the velocity, but even though the streamfunction has fixed Dirichlet boundary conditions, the velocity does not have the same. They do inherhet no-normal flow but that will come automatically from the streamfunction.

Also, thanks for pointing out how to impose boundary conditions more generally, that will come in handy.

I am happy to say that following your help, I do get that part of it running. Now to get this working with the time stepping.

Cheers!

from shenfun.

I am having an odd problem and not sure what to make of it. Below is a samle code that produces the error, albeit in a round about way.

The essence of it is that if I define qk as before to be the div and gradient, then afterwards I cannot compute qk to be the forward transform of qk. I realize this might be a bit odd but I wonder what's going on here? If I redefine qk then it seems to work but this seems less than idea and wonder what I am not understanding here.

`# sample code

a =  1.
P0   = Basis(N[0], family=family, domain=(-L[0]/2,L[0]/2), bc=(-a, a), scaled=True)
P1   = Basis(N[0], family=family, domain=(-L[0]/2,L[0]/2))
F0   = Basis(N[1], family='F', dtype='d', domain=(0., L[1]))
TF0  = TensorProductSpace(comm, (P0, F0), axes=(0, 1))
TF1  = TensorProductSpace(comm, (P1, F0), axes=(0, 1))
TV1  = VectorTensorProductSpace(TF1)

X    = TF0.local_mesh(True)

pr,  pk  = Array(TF0), Function(TF0)                          # Streamfunction
qr,  qk  = Array(TF1), Function(TF1)                          # Potential Vorticity

qr[:] = -2.*np.tanh(X[0])/pow(np.cosh(X[0]),2)        
qk    = TF1.forward(qr, qk)

pr[:] = np.tanh(X[0]) 
pk    = TF0.forward(pr, pk)
qk    = div(grad(pk))                                     # this seems to change qk
qr   = project(qk,TF1).backward(qr)
qk = TF1.forward(qr, qk)

`

from shenfun.

qk is no longer a Function and cannot be used in the last forward call. qk is now an Expr instance. That's just the way regular Python assignments work.

from shenfun.

Thanks. I didn't realize that div(grad) was an expression and not a function. Now I know, and hopefully I won't forget anytime soon.

from shenfun.