<div class="gmail_quote">On Thu, Sep 22, 2011 at 17:21, Milan Mitrovic <span dir="ltr">&lt;<a href="mailto:milan.v.mitrovic@gmail.com">milan.v.mitrovic@gmail.com</a>&gt;</span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div class="im">&gt; Don&#39;t destroy the matrix, just call the preallocation routines when the<br>
&gt; nonzero pattern changes.<br>
&gt; You can actually destroy the matrices from C, but I don&#39;t think there is a<br>
&gt; pure Fortran way to do it (it needs another pointer indirection) and you<br>
&gt; should never need to anyway.<br>
<br>
</div>How should I create a matrix without supplying the preallocation? just<br>
set everything to zero and then call the prealloc routines later?<br></blockquote><div><br></div><div>Just create the matrix and set the sizes, but don&#39;t preallocate. The wrapper you are calling is just a convenience function:</div>
<div><br></div><div><div>PetscErrorCode  MatCreateMPIAIJ(MPI_Comm comm,PetscInt m,PetscInt n,PetscInt M,PetscInt N,PetscInt d_nz,const PetscInt d_nnz[],PetscInt o_nz,const PetscInt o_nnz[],Mat *A)</div><div>{</div><div>  PetscErrorCode ierr;</div>
<div>  PetscMPIInt    size;</div><div><br></div><div>  PetscFunctionBegin;</div><div>  ierr = MatCreate(comm,A);CHKERRQ(ierr);</div><div>  ierr = MatSetSizes(*A,m,n,M,N);CHKERRQ(ierr);</div><div>  ierr = MPI_Comm_size(comm,&amp;size);CHKERRQ(ierr);</div>
<div>  if (size &gt; 1) {</div><div>    ierr = MatSetType(*A,MATMPIAIJ);CHKERRQ(ierr);</div><div>    ierr = MatMPIAIJSetPreallocation(*A,d_nz,d_nnz,o_nz,o_nnz);CHKERRQ(ierr);</div><div>  } else {</div><div>    ierr = MatSetType(*A,MATSEQAIJ);CHKERRQ(ierr);</div>
<div>    ierr = MatSeqAIJSetPreallocation(*A,d_nz,d_nnz);CHKERRQ(ierr);</div><div>  }</div><div>  PetscFunctionReturn(0);</div><div>}</div></div><div><br></div><div><br></div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">

<div class="im"><br>
&gt; What do you mean &quot;complains&quot;? You can certainly use a different<br>
&gt; preconditioning matrix.<br>
<br>
</div>I dont want to use a different matrix... I want to use the same... I<br>
just dont want to get the error that I pasted in the first message...<br>
and I have no idea why I get it in the first place...<br>
</blockquote></div><br><div>A hanging reference was created so the same object was destroyed multiple times. I&#39;m kinda surprised there wasn&#39;t a check earlier to verify that a valid matrix was present. In any case, you shouldn&#39;t create a new matrix inside the Jacobian evaluation function, at least not when calling from Fortran.</div>