Hi Justin, On 15/08/16 02:15, Justin Chang wrote:
Hi all,
So I am implementing PETSc/TAO's semi smooth solver to ensure non-negative solutions for anisotropic diffusion. Attached is the code.
Basically, the call-back functions require a Jacobian and a gradient/constraint evaluation routine. It begins around line 160. Run the code as:
python SS_ex2.py 30 <guess>
where <guess> is either use zeros as initial guess (0) or use the standard formulation solution (i.e., no bounds) as the initial guess (1). The implementation somewhat works if <guess> is 1 or an initial non-zero solution is provided. If no guess is provided, I get a funky solution where everywhere is zero except for the non-homogeneous DirichletBCs.
Also, even if I have an initial guess (1), the values of the BC aren't exactly accurate either (I should have a max value of 1 at the very center of the top face, but I get 0.983941)
Something tells me I am not implementing the BCs correctly. Can anyone tell me what I am possibly doing wrong here?
This is entirely possible. I haven't looked at your formulation in detail, but here's one thing to remember about the way boundary conditions are implemented (summarised from http://firedrakeproject.org/boundary_conditions.html). The particular thing to note is that an operator assembled with boundary conditions will just be the identity on the boundary node rows, because we forward-substitute to the right hand side to make the operator symmetric. If you're doing the linear problem by hand somehow, then you have to do this substitution. You can see this happening in the LinearSolver (in linear_solver.py) for example. Plausibly this means that your gradient should have zero boundary conditions applied to it (do bc.zero(...))? Lawrence