Hi everybody,
I'm building a Poisson solver for a cfd code in a structured grid using fortran
90. Substantially I need to solve a linear system in a cycle with a fixed
matrix and a variable rhs.
I'm working with DAs and (being new to Petsc!) I have a problem.
I have built my matrix with :
...
call DACreate3d(...)
....
call DAGetMatrix(da,MATAIJ,A,err)
......
call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,err)
call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,err)
......
and I've got a matrix A ordered in a Petsc ordering: rows and variables have a
sequential numbering relative to each processor.
On the other hand building the rhs with
......
call DACreateLocalVector(da,qloc,err)
call DAGlobalToLocalBegin(da,q,INSERT_VALUES,qloc,err)
call DAGlobalToLocalEnd(da,q,INSERT_VALUES,qloc,err)
call VecGetArrayF90(qloc,qloc_a,err)
cont=0
do k=gzs,gzs+gzm-1
do j=gys,gys+gym-1
do i=gxs,gxs+gxm-1
if ( ...... ) then
cont=cont+1
qloc_a(cont)=....
else
......
endif
enddo
enddo
enddo
call VecRestoreArrayF90(qloc,qloc_a,err)
call DALocalToGlobal(da,qloc,INSERT_VALUES,q,err)
I end up with a vector q ordered in the natural ordering (as if it is built
with one processor). The questions are:
Is there a way to end up with the same ordering for both the matrix and the rhs
?
If the only way is to use the AO mapping obtained with DAGetAO(), is there an
example code that shows how to use the resulting ao to remap the vector (or the
matrix) in the Petsc (or the natural) order?
I've seen that I can have the remapping indices but then I haven't understood
how to use them to remap the vector (or the matrix).
Regards
Valerio Grazioso
-------------- next part --------------
An HTML attachment was scrubbed...
URL:
<http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20100622/7f11067f/attachment-0001.htm>
