<html><head><meta http-equiv="Content-Type" content="text/html; charset=utf-8"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class="">Hello Stefano,<br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On 20 Oct 2018, at 9:22 AM, Stefano Zampini <<a href="mailto:stefano.zampini@gmail.com" class="">stefano.zampini@gmail.com</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="auto" class="">Pierre,<div dir="auto" class=""><br class=""></div><div dir="auto" class="">Thanks to report this. I think this should be fixed in MatHermitianTransposeAdd. Would you like to open a PR?</div></div></div></blockquote><div><br class=""></div><div>Done.</div><div>Thanks,</div><div>Pierre</div><br class=""><blockquote type="cite" class=""><div class=""><div dir="auto" class=""><div dir="auto" class="">Thanks</div><div dir="auto" class="">Stefani</div></div><br class=""><div class="gmail_quote"><div dir="ltr" class="">Il Sab 20 Ott 2018, 08:16 Pierre Jolivet <<a href="mailto:pierre.jolivet@enseeiht.fr" class="">pierre.jolivet@enseeiht.fr</a>> ha scritto:<br class=""></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hello,<br class="">
I’m using a MatNest with some blocks made of matrices created with MatCreateHermitianTranspose. I end up with the following trace when doing a simple MatMult.<br class="">
Should this be fixed in matnest.c, matrix.c or htransm.c?<br class="">
<br class="">
Thanks,<br class="">
Pierre<br class="">
<br class="">
(lldb) down<br class="">
frame #11: 0x000000011ad06830 libpetsc.3.010.dylib`MatMult_Nest(A=0x00007fcc2eb50260, x=0x00007fcc2f019460, y=0x00007fcc2f01ca60) at matnest.c:52<br class="">
49 for (j=0; j<nc; j++) {<br class="">
50 if (!bA->m[i][j]) continue;<br class="">
51 /* y[i] <- y[i] + A[i][j] * x[j] */<br class="">
-> 52 ierr = MatMultAdd(bA->m[i][j],bx[j],by[i],by[i]);CHKERRQ(ierr);<br class="">
53 }<br class="">
54 }<br class="">
55 for (i=0; i<nr; i++) {ierr = VecRestoreSubVector(y,bA->isglobal.row[i],&by[i]);CHKERRQ(ierr);}<br class="">
(lldb) <br class="">
frame #10: 0x000000011a9d6399 libpetsc.3.010.dylib`MatMultAdd(mat=0x00007fcc2ebfde60, v1=0x00007fcc2f029c60, v2=0x00007fcc2f024660, v3=0x00007fcc2f024660) at matrix.c:2517<br class="">
2514 if (!mat->ops->multadd) SETERRQ1(PetscObjectComm((PetscObject)mat),PETSC_ERR_SUP,"No MatMultAdd() for matrix type '%s'",((PetscObject)mat)->type_name);<br class="">
2515 ierr = PetscLogEventBegin(MAT_MultAdd,mat,v1,v2,v3);CHKERRQ(ierr);<br class="">
2516 ierr = VecLockPush(v1);CHKERRQ(ierr);<br class="">
-> 2517 ierr = (*mat->ops->multadd)(mat,v1,v2,v3);CHKERRQ(ierr);<br class="">
2518 ierr = VecLockPop(v1);CHKERRQ(ierr);<br class="">
2519 ierr = PetscLogEventEnd(MAT_MultAdd,mat,v1,v2,v3);CHKERRQ(ierr);<br class="">
2520 ierr = PetscObjectStateIncrease((PetscObject)v3);CHKERRQ(ierr);<br class="">
(lldb) <br class="">
frame #9: 0x000000011ad136d0 libpetsc.3.010.dylib`MatMultAdd_HT(N=0x00007fcc2ebfde60, v1=0x00007fcc2f029c60, v2=0x00007fcc2f024660, v3=0x00007fcc2f024660) at htransm.c:24<br class="">
21 PetscErrorCode ierr;<br class="">
22 <br class="">
23 PetscFunctionBegin;<br class="">
-> 24 ierr = MatMultHermitianTransposeAdd(Na->A,v1,v2,v3);CHKERRQ(ierr);<br class="">
25 PetscFunctionReturn(0);<br class="">
26 }<br class="">
27 <br class="">
(lldb) <br class="">
frame #8: 0x000000011a9d88e5 libpetsc.3.010.dylib`MatMultHermitianTransposeAdd(mat=0x00007fcc2e977260, v1=0x00007fcc2f029c60, v2=0x00007fcc2f024660, v3=0x00007fcc2f024660) at matrix.c:2629<br class="">
2626 ierr = MatMultTranspose(mat,w,z);CHKERRQ(ierr);<br class="">
2627 ierr = VecDestroy(&w);CHKERRQ(ierr);<br class="">
2628 ierr = VecConjugate(z);CHKERRQ(ierr);<br class="">
-> 2629 ierr = VecWAXPY(v3,1.0,v2,z);CHKERRQ(ierr);<br class="">
2630 ierr = VecDestroy(&z);CHKERRQ(ierr);<br class="">
2631 }<br class="">
2632 ierr = VecLockPop(v1);CHKERRQ(ierr);<br class="">
(lldb) <br class="">
frame #7: 0x000000011a8f33b5 libpetsc.3.010.dylib`VecWAXPY(w=0x00007fcc2f024660, alpha=1, x=0x00007fcc2f024660, y=0x00007fcc2f039260) at rvector.c:795<br class="">
792 VecCheckSameSize(x,3,y,4);<br class="">
793 VecCheckSameSize(x,3,w,1);<br class="">
794 if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");<br class="">
-> 795 if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");<br class="">
796 PetscValidLogicalCollectiveScalar(y,alpha,2);<br class="">
797 <br class="">
798 ierr = PetscLogEventBegin(VEC_WAXPY,x,y,w,0);CHKERRQ(ierr);<br class="">
<br class="">
</blockquote></div>
</div></blockquote></div><br class=""></body></html>