On Thu, Oct 10, 2013 at 2:33 PM, Jed Brown <[email protected]> wrote:
> Christophe Ortiz <[email protected]> writes: > > I am facing some problems with this when filling in the Jacobian for a > > problem of 2 diffusion equations (dof=2). See ex25.c. > > > > If I'm right, then for node i, there should be 3 blocks of 2x2: > > > > X 0 | Y 0 | X 0 > > 0 X | 0 Y | 0 X > > > > Is it correct ? > > > > First I tried what you said, to pass a one dimensional array: > > row = i; > > cols[0] = i-1; cols[1] = i; cols[2] = i+1; > > > > j = 0; > > vallin[j] = X; vallin[j++] = 0; > > The expression j++ in vallin[j++] is post-increment, so the two > statements above both write to vallin[0] and vallin[11] is undefined > (whatever garbage was in the array before). > Thanks. I just found out few minutes ago :-). Hehe. The point is that I come from 15 years of Fortran...I have to refresh the C that I learnt at Faculty and that is deeply buried in my memory ! > > > Then, I tried a 2-dimensional array. I saw in some examples of > ts/tutorials > > that it is also possible to pass a multi-dimensional arrays. After all, > in > > memory the rows are stored contiguously, so I guess it is equivalent to a > > one-dimensional array. > > > > j = 0; k = 0; > > valdiff[j][k] = X; valdiff[j][k++] = 0; > > Same problem here. > > > valdiff[j][k++] = Y; valdiff[j][k++] = 0; > > valdiff[j][k++] = X; valdiff[j][k++] = 0; > > > > j= 1; k = 0; > > valdiff[j][k] = 0; valdiff[j][k++] = X; > > valdiff[j][k++] = 0; valdiff[j][k++] = Y; > > valdiff[j][k++] = 0; valdiff[j][k++] = X; > > > > ierr = > > MatSetValuesBlocked(*Jpre,1,&row,3,cols,&valdiff[0][0],INSERT_VALUES); > > > > I get the same error message. > > What did I do wrong ? > > > > In ex25.c they use a 3-dimensional array: > > > > const PetscScalar vals[2][3][2] = {{{dxxL,0},{a *hx+dxx0,0},{dxxR,0}}, > > > > {{0,dxxL},{0,a*hx+dxx0},{0,dxxR}}}; > > > > and &vals[0][0][0] is passed in MatSetValuesBlocked (....,...., > > &vals[0][0][0],...), and it works, though I don't understand how. > > That grouping in the 3D array is the same ordering, but (perhaps) makes > it easier to see which expression refers to which entry of the matrix. > Yes, probably.
