Thanks Stefano, I will try what you suggest. Matt - my DM is a composite between the redundant field (loading coefficient, which is included in the Newton solve in Riks' method) and the displacements, which are represented by a 3D DA with 3 dof. I am using finite difference.
Probably my problem comes from confusion over how the composite DM is organized. I am using FormFunction(), and within that I call DMCompositeGetLocalVectors(), DMCompositeScatter(), DMDAVecGetArray(), and for the Jacobian, DMCompositeGetLocalISs() and MatGetLocalSubmatrix() to split J into Jbb, Jbh, Jhb, and Jhh, where b is the loading coefficient, and h is the displacements). The values of each submatrix are set using MatSetValuesLocal(). I'm most suspicious of the part of the Jacobian routine where I calculate the rows of Jhb, the columns of Jbh, and the corresponding values. I take the DA coordinates and ix,iy,iz, then calculate the row of Jhb as ((((iz-info->gzs)*info->gym + (iy-info->gys))*info->gxm + (ix-info->gxs))*info->dof+c), where info is the DA local info and c is the degree of freedom. The same calculation is performed for the column of Jbh. I suspect that the indexing of the DA vector is not so simple, but I don't know for a fact that I'm doing this incorrectly nor how to do this properly. Thanks for all the help! On Nov 9, 2017 8:44 AM, "Matthew Knepley" <[email protected]> wrote: On Thu, Nov 9, 2017 at 12:14 AM, zakaryah . <[email protected]> wrote: > Well the saga of my problem continues. As I described previously in an > epic thread, I'm using the SNES to solve problems involving an elastic > material on a rectangular grid, subjected to external forces. In any case, > I'm occasionally getting poor convergence using Newton's method with line > search. In troubleshooting by visualizing the residual, I saw that in data > sets which had good convergence, the residual was nevertheless > significantly larger along the boundary between different processors. > Likewise, in data sets with poor convergence, the residual became very > large on the boundary between different processors. The residual is not > significantly larger on the physical boundary, i.e. the global boundary. > When I run on a single process, convergence seems to be good on all data > sets. > > Any clues to fix this? > It sounds like something is wrong with communication across domains: - If this is FEM, it sounds like you are not adding contributions from the other domain to shared vertices/edges/faces - If this is FDM/FVM, maybe the ghosts are not updated What DM are you using? Are you using the Local assembly functions (FormFunctionLocal), or just FormFunction()? Thanks, Matt -- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.caam.rice.edu/~mk51/>
