That&#39;s great. Thank you very much, Jed :).<br><br>Best,<br>Yujie<br><br><div class="gmail_quote">On Sat, Feb 11, 2012 at 11:06 AM, Jed Brown <span dir="ltr">&lt;<a href="mailto:jedbrown@mcs.anl.gov">jedbrown@mcs.anl.gov</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"><div class="gmail_quote">On Sat, Feb 11, 2012 at 10:58, recrusader <span dir="ltr">&lt;<a href="mailto:recrusader@gmail.com" target="_blank">recrusader@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>    &quot;    ierr = MatCreateMPIAIJ (libMesh::COMM_WORLD,<br>                              </div><div><div>  m_local, n_local,<br>
                                m_global, n_global,<br>
                                PETSC_NULL, (int*) &amp;n_nz[0],<br>                                PETSC_NULL, (int*) &amp;n_oz[0], &amp;_mat);<br>             CHKERRABORT(libMesh::COMM_WORLD,ierr);<br><br>      MatSetOption(_mat,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE); //by Yujie<br>



      std::cout&lt;&lt;&quot;MatSetOption&quot;&lt;&lt;std::endl;&quot;<br><br>I 
run the same codes in CPU and GPU modes (the same parameters except that
 GPU uses &#39;-vec_type mpicusp -mat_type mpiaijcusp&#39;). I can find 
&quot;MatSetOption&quot; output from both the modes. Does that mean that the codes
 set the options for both the modes?</div></div></blockquote></div><div><br></div></div>Libmesh is calling MatSetFromOptions() after MatCreateMPIAIJ() which means the preallocation information will be lost if the type is changed. The code should be written differently<div>

<br></div><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>  ierr = MatSetType(*A,MATAIJ);CHKERRQ(ierr);</div>

<div>  ierr = MatSetOptionsPrefix(*A,optional_prefix);CHKERRQ(ierr);</div><div>  ierr = MatSetFromOptions(*A);CHKERRQ(ierr);</div><div>  ierr = MatMPIAIJSetPreallocation(*A,d_nz,d_nnz,o_nz,o_nnz);CHKERRQ(ierr);</div><div>

  ierr = MatSeqAIJSetPreallocation(*A,d_nz,d_nnz);CHKERRQ(ierr);</div><div><br></div></div><div>I can talk to the libmesh developers about making this change.</div>
</blockquote></div><br>