-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 Hi Justin, On 15/09/15 09:49, Justin Chang wrote:
Lawrence,
Thanks for all of your help. Here’s what I ended up with: —————— def ObjGrad_D(tao, petsc_x, petsc_g): with b_D_nobc.dat.vec as b_D_vec: A_D_nobc.M.handle.mult(petsc_x,petsc_g) xtHx = petsc_x.dot(petsc_g) xtf = petsc_x.dot(b_D_vec) petsc_g.axpy(-1.0,b_D_vec) arr = petsc_g.array for bc in bc_D: nodes = bc.nodes[bc.nodes < bc.function_space().dof_dset.size] if isinstance(bc.function_arg, Function): with bc.function_arg.dat.vec_ro as bc_vec: arr[nodes] = bc_vec.array_r[nodes] else: arr[nodes] = float(bc.function_arg) return 0.5*xtHx - xtf ——————— However it seems something is still off about the boundary conditions I am applying. I have attached three different scenarios for my two solutions psi_A and psi_B. The boundary conditions are:
psi_A: set to 1 at x == 0.0 and 0.0 <= y <= 0.5, 0 everywhere else psi_B: set to 1 at x == 0.0 and 0.5 < y <= 1.0, 0 everywhere else
- no_opt_psiA/B.png illustrates solving the standard diffusion without any optimization. As far as I know, these solutions are “correct” in terms of the boundary conditions.
- opt_wo_bounds_psiA/B.png illustrates solving the standard diffusion with TAO’s optimization but without any variable bounds. I would expect this to have the same solution as the no_opt cases but it seems to me the boundary conditions are messed up; face ID’s 3 and 4 (i.e., top and bottom) should have homogeneous bcs but it seems they are not registering.
- opt_w_bounds_psiA/B.png now includes the variable bounds set by these lines:
lb = Function(W) ub = Function(W) lb.assign(0) ub.assign(1) tao = PETSc.TAO().create(MPI.COMM_WORLD) with lb.dat.vec as lb_vec, ub.dat.vec as ub_vec: tao.setVariableBounds(lb_vec,ub_vec)
The boundary conditions are also screwed up at the top and bottom faces. In this case, psi_B (corresponding to W.sub(1)) is even more screwed up - the non-homogeneous part of psi_B isn’t even registering (hence concentrations of < 1).
Replacing the RHS of the arr[nodes] assignments with 0 didn’t change the solution either. Any idea what’s going on here?
I don't have any immediate insights, sorry. But some comments. When you're minimising 0.5 x^T H x - x^T f What parts of that do you expect to obey the boundary conditions? I presume that f should? I'm a bit lost in all the code. Can you maybe write down the maths (with the specification of the bcs) and we can go from them instead? Lawrence -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 iQEcBAEBAgAGBQJV9+zQAAoJECOc1kQ8PEYvex8IAI+biiyYkJmhE57nyS13tiNs PxS/TLu1UFMAHlXs272S3r+U4DMXeVFHeb2kAf9hwSUf9g25+38rAfSYX2SH5KX2 +uXpCeOJAuxquS/zI9i7y8AP+ftbTMbmJOGLjnFZRlEDsHVjvqEToPzO13JnXCYp SIpGId6aO0FbFwT94ZDIV9fQKn2bYsLw9FSjWAEVCXoZiL6yEeyNstH32jSf8EU4 JrXsyekL4IGc8t5UtfpWiumVRMwRdfQntA+oqpEyeKiqDLFdTACAMzuAxDyXwbXH GFLcN3yYMU52R82KpdsymOEodYRaliF4CvULXGmEQwJkmoXJJvDUFnPF0b37VJ0= =JH/J -----END PGP SIGNATURE-----