[petsc-users] DA matrices and vectors ordering

[Valerio Grazioso]

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>