Specifying boundary conditions for both fields in a coupled PDE

[davenpi]

Hi all,

This package seems great! I've been looking for something like this for a while so I'm happy to try it out.

I want to simulate the PDE pictured below.

The sigma field is just some source term that I specify manually. (added equation for sigma so PDE class would recognize variable)

When I integrate this equation

I find that the left boundary condition on the field u is not satisfied (didn't check derivative BC on right side yet). I want u(0) = 0 but when I inspect the solution fields I find that u(0) is not actually zero.

It seems that the BCs have been ignored or I haven't specified them properly. I tried making a tuple for each field
bc = [({"value":0},{"derivative": -0.047}), ({"value":0},{"value":10})] (first tuple has u BCs and second tuple has sigma BCs) but when I do that I get the following error

I've looked through the docs for info on specifying BCs for both fields in a coupled PDE but I wasn't able to find any. Do you know how I can solve this issue?

[david-zwicker]

Thank you for the report. I think you specified the boundary conditions correctly. In fact, the image you posted looks as if u(0)=0 is satisfied. Note that the first support point is at x=0.025 in your case, because the cell centers are used to store field values; see documentation. If you want to obtain the value at x=0, you can use sol.interpolate(0).

[davenpi]

Thanks for the response. I didn't realize the grid was using the center point of intervals. I am confused about what specifying
bc = [({"value": 0}, {"derivative": -0.047})] means. Taking the center point of the first interval as the start, my understanding was that this would mean u(x=0.025) = 0. But here we don't have that condition satisfied. Then when I look at the interpolation at x=0 it gives the same value as at x=0.025

The docs say the point argument in the interpolate method is given in grid coordinates. It seems that means sol[0].interpolate(0) computes the interpolation at x=0 so why does sol[0].interpolate(0) give the same value as sol[0](x=0.025) ?

And finally, how would I specify different boundary conditions for each field? I couldn't find any description of how to do that in the docs.

Thanks again and let me know if I should clarify anything.

[david-zwicker]

I converted the issue to a discussion, since this format lends itself better for the discussion we're having and it is also archived. Feel free to open an independent issue if you find a specific problem with the package.

While the grid points are specified at the cell centers, the boundary conditions are still imposed at the boundary, which explains why the value is not zero at the first grid point. I now noticed that the interpolation at the boundary does not extrapolate the field, but uses the value of the closest grid point, which is a good general default. However, it would be nice if the interpolator would know about the boundary conditions, which would allow it to determine the value at the boundary correctly. This feature might come in the future, but I currently do not have any time to implement it.

Considering the boundary conditions for different fields, you need to use the bc_ops argument of the PDE class. In your case, you could use

pde.PDE(
    {'u': '...', 'sig_field': '...'},
    bc_ops={'u:*': ({'value': 0, 'derivative': 1), 'sig_field:*': {'derivative': 0}}
)

However, an even better solution might be to not introduce the sig_field and instead use the consts arguments of pde.PDE to supply the gradient of this field directly, so it does not need to be calculated each time step (and the field is not advanced in time).

[davenpi]

Okay thanks for the info. Happy to move this to a discussion. I was expecting a linear fit since the interpolate method parameter has a default value "linear" but, as you point out, zero slope is linear! It seems setting the bc_ops as you suggested ended up working too.

[david-zwicker]

That's good to hear! I added support for interpolation with more flexible boundary conditions to the TODO list, but I have no idea when I'll find the time to implement this.