Updating BC Expression

variational formulation
1

Hello,

I want to update my boudnary condition and I’m running into an error which, I think is due to a worng definition or expresison, I’m not sure. However, I’m not sure how to solve this problem.

pres_oo = Expression('art_pressure', art_pressure=0., degree=1)

# Define boundary conditions
bcu_noslip1  = DirichletBC(V, Constant((0,0)), walls)
bcp_inflow  = DirichletBC(Q, Constant(8), inflow)
bcp_outflow = DirichletBC(Q, pres_oo, outflow)

…

Rp = 1e-6
C = .5 * Rp
Z = (1.1 * Rp) * 0.05
RC = Rp * C 
Q_outlet=0.

P_out_new = Expression('((Z+(Rp*dt)/(RC+dt)-((Z*RC)/(RC+dt))*Q_outlet+((RC)/(RC+dt))*p_n))',Z=Z,Rp=Rp,RC=RC,dt=dt,Q_outlet=Q_outlet,p_n=p_n,degree=1)

…

pres_oo.art_pressure=P_out_new #update BC
This update BC part is in the time-stepping part.
The above snippets are the relevant parts of my code. The idea is that every time step the new P_out is calculated and applied as BC. However, when running this code i get:

Traceback (most recent call last):
  File "tut09.py", line 194, in <module>
    pres_oo.art_pressure=P_out_new      #update BC
  File "/usr/local/lib/python3.6/dist-packages/dolfin/function/expression.py", line 439, in __setattr__
    self._parameters[name] = value
  File "/usr/local/lib/python3.6/dist-packages/dolfin/function/expression.py", line 299, in __setitem__
    self._cpp_object.set_property(key, value)
RuntimeError: No such property

I cannot seem to figure out whats going wrong, when I want to print the value of P_out_new during the time step it does not give me a value but f_55 or f_54.

Thanks in advance!
Cheers,
Jack

2

Hi @jacktat
Please consider using UserExpression, as illustrated below:

from dolfin import *
mesh = UnitSquareMesh(10,10)
V = FunctionSpace(mesh, "CG", 1)

class press_oo(UserExpression):
    def __init__(self, art_pressure, **kwargs):
        super().__init__(kwargs)
        self.art_pressure = art_pressure
    def eval(self, values, x):
        values[0] = self.art_pressure(x)
    def value_shape(self):
        return ()

art_pressure = Function(V)
pres_oo = press_oo(art_pressure)
print(assemble(pres_oo*dx(domain=mesh)))
x = SpatialCoordinate(mesh)
pres_oo.art_pressure.assign(project(x[0]*x[1],V))

print(assemble(pres_oo*dx(domain=mesh)))

which yields the output

0.0
0.24999999999999983
3

Thanks for the response Dokken. How does a userExpression look like for the P_out_new? I figured this also has to become a userExpression where I update Q_outlet every timestep. But I don’t know how how to get that full formula
P_out_new = Expression('((Z+(Rp*dt)/(RC+dt)-((Z*RC)/(RC+dt))*Q_outlet+((RC)/(RC+dt))*p_n))',Z=Z,Rp=Rp,RC=RC,dt=dt,Q_outlet=Q_outlet,p_n=p_n,degree=1)

as a userExpression.

any tips?

4

It is quite straight-forward:

from dolfin import *
mesh = UnitSquareMesh(10,10)
V = FunctionSpace(mesh, "CG", 1)

Rp = 1e-6
C = .5 * Rp
Z = (1.1 * Rp) * 0.05
RC = Rp * C
dt = 0.1

class p_out_new(UserExpression):
    def __init__(self, p_n, q_outlet,**kwargs):
        super().__init__(kwargs)
        self.p_n = p_n
        self.q_outlet = q_outlet
    def eval(self, values, x):
        values[0] = (Z + (Rp * dt) / (RC + dt)
                     - (Z * RC)/(RC + dt)*self.q_outlet(x)
                     + RC / (RC + dt) * self.p_n(x))

    def value_shape(self):
        return ()

p_new = p_out_new(Function(V), Function(V))
print(assemble(p_new*dx(domain=mesh)))
x = SpatialCoordinate(mesh)



p_new.p_n.assign(project(x[0]*x[1],V))
p_new.q_outlet.assign(project(x[1],V))
print(assemble(p_new*dx(domain=mesh)))
1 Like