<div dir="ltr">Hi,<div><br></div><div>I am having trouble calling KSPSolveTranspose under MATSOLVERPASTIX.</div><div>The following is a very simple C code solving A^T x = b to illustrate the issue. The matrix A and vector b are built following the PETSC tutorial.</div>
<div>KSPSolve works fine but when it comes to KSPSolveTranspose, it spits out the error:</div><div><br></div><div><div><font face="courier new, monospace">[0]PETSC ERROR: --------------------- Error Message ------------------------------------</font></div>
<div><font face="courier new, monospace">[0]PETSC ERROR: No support for this operation for this object type!</font></div><div><font face="courier new, monospace">[0]PETSC ERROR: Matrix type mpiaij!</font></div></div><div>
<br></div><div>The code runs fine if I change the solver package <span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MATSOLVERPASTIX </span><span style="color:rgb(0,0,0);font-size:14px"><font face="arial, helvetica, sans-serif">to </font><font face="courier new, monospace">MATSOLVERMUMPS</font><font face="arial, helvetica, sans-serif">.</font></span></div>
<div><span style="color:rgb(0,0,0);font-size:14px"><font face="arial, helvetica, sans-serif">However, I need the PASTIX solver and the transpose solution for solving a large problem.</font></span></div><div><span style="color:rgb(0,0,0);font-size:14px"><font face="arial, helvetica, sans-serif">Could you please let me know what is wrong with the PASTIX solver or with my usage?</font></span></div>
<div><span style="color:rgb(0,0,0);font-size:14px"><font face="arial, helvetica, sans-serif"><br></font></span></div><div><span style="color:rgb(0,0,0);font-size:14px"><font face="arial, helvetica, sans-serif">Best regards,</font></span></div>
<div><span style="color:rgb(0,0,0);font-size:14px"><font face="arial, helvetica, sans-serif">Zin</font></span></div><div><span style="color:rgb(0,0,0);font-size:14px"><font face="arial, helvetica, sans-serif"><br></font></span></div>
<div><br></div><div><font face="courier new, monospace">/******************/</font></div><div><div><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">#include <stdlib.h></span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">#include <petsc.h></span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">#include <string.h></span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">#include <complex.h></span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">int main(int argc, char **argv)</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">{</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> /* -------Initialize ------*/</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">PetscInitialize(&argc, &argv, PETSC_NULL, PETSC_NULL);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">PetscPrintf(PETSC_COMM_WORLD,"--------Initializing------ \n");</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">PetscErrorCode ierr;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">Vec x,y;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">Mat A;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">int N=10;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">int low, high, i;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">double val;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">int start, end;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">double v[3];</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">int row, cols[3];</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">ierr = VecCreateMPI(PETSC_COMM_WORLD, PETSC_DECIDE,N, &x);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">ierr = VecDuplicate(x,&y); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">VecGetOwnershipRange(x, &low, &high);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">val = low*10.0;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">for(i = low; i < high; ++i) {</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> VecSetValues(x, 1, &i, &val, INSERT_VALUES);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> val += 10.0;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">}</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">VecAssemblyBegin(x);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">VecAssemblyEnd(x);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MatCreate(PETSC_COMM_WORLD, &A);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MatSetType(A,MATMPIAIJ);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MatSetSizes(A,PETSC_DECIDE, PETSC_DECIDE, N, N);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MatMPIAIJSetPreallocation(A, 3, PETSC_NULL, 3, PETSC_NULL);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">v[0] = -1.0; v[1] = 2.0; v[2] = -1.0;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MatGetOwnershipRange(A,&start,&end);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">for(row = start; row < end; row++) {</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> cols[0] = row-1; cols[1] = row; cols[2] = row+1;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">if (row == 0) {</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> MatSetValues(A,1,&row,2,&cols[1],&v[1],INSERT_VALUES);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">} else if (row == N-1) {</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> MatSetValues(A,1,&row,2,cols,v,INSERT_VALUES);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">} else {</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> MatSetValues(A,1,&row,3,cols,v,INSERT_VALUES);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">}</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">}</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">PetscPrintf(PETSC_COMM_WORLD,"A and x setup. \n");</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">KSP ksp;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">PC pc; </span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPSetType(ksp, KSPGMRES);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = PCSetType(pc,PCLU);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = PCFactorSetMatSolverPackage(pc,MATSOLVERPASTIX);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = PCSetFromOptions(pc);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPSetTolerances(ksp,1e-14,PETSC_DEFAULT,PETSC_DEFAULT,20);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPSetOperators(ksp,A,A,SAME_NONZERO_PATTERN);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPSolve(ksp,x,y);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">double p;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = VecDot(x,y,&p); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">PetscPrintf(PETSC_COMM_WORLD,"x dot y: %g. \n", p);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">Vec u;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">ierr = VecDuplicate(x,&u); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">ierr = KSPSolveTranspose(ksp,y,u);CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">PetscPrintf(PETSC_COMM_WORLD,"Transpose solved. \n");</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = VecDestroy(&x); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = VecDestroy(&y); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = MatDestroy(&A); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = KSPDestroy(&ksp); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> {</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">   int rank;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">   MPI_Barrier(PETSC_COMM_WORLD);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> }</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> </span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"> ierr = PetscFinalize(); CHKERRQ(ierr);</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">return 0;</span><br style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">
<span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px">}</span><br></div><div><span style="font-family:'courier new',monospace">/******************/</span><span style="color:rgb(0,0,0);font-family:'Courier New',Courier,monospace;font-size:14px"><br>
</span></div>-- <br><div dir="ltr">Zin Lin<br><br></div>
</div></div>