<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Mon, Sep 18, 2017 at 6:16 PM, David Gross <span dir="ltr"><<a href="mailto:davegwebb10@gmail.com" target="_blank">davegwebb10@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><div><div>Hi Matt,<br></div>I took a deeper look into the dev documentation and the src/mat code. </div></div></div></div></blockquote><div><br></div><div>All the source is online, so I think that is a good way to talk about it in email.</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><div>From what I understand there are multiple places that the code needs to be added or referenced in. I created a copy of AXPY in axpy.c </div></div></div></div></blockquote><div><br></div><div>Lets look at MatAXPY():</div><div><br></div><div> <a href="https://bitbucket.org/petsc/petsc/src/a2d7f5eaef7c3b636ab943ac50d4f7c86905df62/src/mat/utils/axpy.c?at=master&fileviewer=file-view-default#axpy.c-22">https://bitbucket.org/petsc/petsc/src/a2d7f5eaef7c3b636ab943ac50d4f7c86905df62/src/mat/utils/axpy.c?at=master&fileviewer=file-view-default#axpy.c-22</a></div><div><br></div><div>You just need to make a copy of this function. You want Z = X * Y, right?, (or maybe Y = X * Y) so its</div><div><br></div><div> MatPointwiseMult(Mat Z, Mat X, Mat Y)</div><div><br></div><div>After the initial checks for compatibility, it dispatches</div><div><br></div><div> <a href="https://bitbucket.org/petsc/petsc/src/a2d7f5eaef7c3b636ab943ac50d4f7c86905df62/src/mat/utils/axpy.c?at=master&fileviewer=file-view-default#axpy.c-36">https://bitbucket.org/petsc/petsc/src/a2d7f5eaef7c3b636ab943ac50d4f7c86905df62/src/mat/utils/axpy.c?at=master&fileviewer=file-view-default#axpy.c-36</a></div><div><br></div><div>Here lets ignore the version optimized for a particular storage format, and just use the generic version</div><div><br></div><div> <a href="https://bitbucket.org/petsc/petsc/src/a2d7f5eaef7c3b636ab943ac50d4f7c86905df62/src/mat/utils/axpy.c?at=master&fileviewer=file-view-default#axpy.c-58">https://bitbucket.org/petsc/petsc/src/a2d7f5eaef7c3b636ab943ac50d4f7c86905df62/src/mat/utils/axpy.c?at=master&fileviewer=file-view-default#axpy.c-58</a></div><div><br></div><div>Here we take advantage of the built-in ADD_VALUES, but there is no MULT_VALUES option. Thus, we would need to</div><div>get the row of X and Y, compare columns (which come out sorted) and then insert the result. However, the alternative is</div><div>have it error if MatStructure is DIFFERENT_NONZERO_PATTERN, so you can assume that both rows are identical and</div><div>dispense with the column checking.</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><div>and compiled without issue and could even reference it in my code without compiler error, but found that it wasn't actually doing anything as best I can tell. Am I correct in that one must make implementations for each matrix type in mat/impls? Also I assume that it needs to be added to matipl.h, include/petscmat.h, and finclude/petscmat.h .<br><br></div>Is there somewhere that documents the process of adding functions to PETSc classes? I couldn't find anything in the developers manual.<br></div></div></div></blockquote><div><br></div><div>There is a discussion of this in the tutorials, and in here <a href="http://www.mcs.anl.gov/petsc/developers/index.html">http://www.mcs.anl.gov/petsc/developers/index.html</a> there is a discussion of adding different types of things.</div><div><br></div><div> THanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div></div>Regards,<br></div>Dave<br></div><div class="gmail_extra"><br><div class="gmail_quote">On Thu, Sep 14, 2017 at 10:07 PM, David Gross <span dir="ltr"><<a href="mailto:davegwebb10@gmail.com" target="_blank">davegwebb10@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><div><div><div>Hi Matt,<br></div>Noted for the naming convention.<br><br></div>Yes, it is a Newton method (see Pan, V. and Reif, J., Fast and Efficient Parallel Solution of Dense Linear Systems, Computers Math. Applications. Vol. 17, No. 11, pp. 1481-1491, 1989)<br><br></div>The dense matrix I have is repeatedly inverted while slowly changing such that the previous inverse is a near perfect guess for the new inverse.<br><br></div>Regards,<br></div>Dave<br></div><div class="gmail-m_-4173110423624862420HOEnZb"><div class="gmail-m_-4173110423624862420h5"><div class="gmail_extra"><br><div class="gmail_quote">On Thu, Sep 14, 2017 at 2:49 AM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><span>On Wed, Sep 13, 2017 at 6:14 PM, David Gross <span dir="ltr"><<a href="mailto:davegwebb10@gmail.com" target="_blank">davegwebb10@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><div><div>Hi Matt,<br></div>Thank you for getting back to me. Your answer confirms what I thought in terms of existing functionality. I think it shouldn't be too hard to make a copy of MatAXPY to MatAXY where it performs Xij = A*Xij*Yij (or without the A). I could then do the MatNorm of the resulting matrix to get what I need.</div><div><br></div><div>Is a MatAXY function desirable as a source contribution?</div></div></div></div></blockquote><div><br></div></span><div>Yes. I would prefer calling it MatPointwiseMult, since you can see it as VecPointwiseMult on a Vec obtained</div><div>from forgetting the linear operator structure of the matrix (forgetful functor).</div><span><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div>I am hoping to use PETSc for performing basic vector and matrix operations on dense matrices and 1D vectors. The main uses are matmult, matmatmult and matrix additions and scaling. The application is for implementing a parallel version on an existing Pan-Reif matrix inversion algorithm.</div></div></div></blockquote><div><br></div></span><div>Is this Newton's method on the vector space of matrices?</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div><div class="gmail-m_-4173110423624862420m_697185988798551954h5"><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div> The choice of using PETSc is mostly due to us already using it in the same program to solve sparse matrices (with MUMPS) with the goal of avoiding adding yet another package (ex ScaLAPACK/PBLAS) into the code even if other packages may be more directly oriented towards my application.<br><br></div>Regards,<br></div>Dave<br><div><div><div><div><div><br><br></div></div></div></div></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Mon, Sep 11, 2017 at 2:19 AM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><span>On Sun, Sep 10, 2017 at 5:51 PM, David Gross <span dir="ltr"><<a href="mailto:davegwebb10@gmail.com" target="_blank">davegwebb10@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><div><div>Hello,<br></div>I was wondering if there was a matrix equivalent to the vecDot function (Frobenius inner product)? As far as I can tell the closest thing is MatNorm with NORM_FROBENIUS, but obviously this is acting on only one matrix.<br><br></div>If there is not a built in function, what is the best way to compute this? I am working fortran90.<br></div></div></div></blockquote><div><br></div></span><div>We do not have this. However, it would be trivial to add since we have</div><div><br></div><div> <a href="http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatAXPY.html" target="_blank">http://www.mcs.anl.gov/petsc<wbr>/petsc-current/docs/manualpage<wbr>s/Mat/MatAXPY.html</a></div><div><br></div><div>since you just replace + with * in our code. You could argue that we should have written for</div><div>a general ring, but C makes this cumbersome. Do you think you could make the change?</div><div><br></div><div>What are you using this for?</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div></div>Regards,<br></div>Dave<span class="gmail-m_-4173110423624862420m_697185988798551954m_-4157210347045451781m_674139435811179689HOEnZb"><font color="#888888"><br></font></span></div><span class="gmail-m_-4173110423624862420m_697185988798551954m_-4157210347045451781m_674139435811179689HOEnZb"><font color="#888888">
</font></span></blockquote></div><span class="gmail-m_-4173110423624862420m_697185988798551954m_-4157210347045451781m_674139435811179689HOEnZb"><font color="#888888"><br><br clear="all"><span class="gmail-HOEnZb"><font color="#888888"><span class="gmail-m_-4173110423624862420m_697185988798551954m_-4157210347045451781HOEnZb"><font color="#888888"><div><br></div>-- <br><div class="gmail-m_-4173110423624862420m_697185988798551954m_-4157210347045451781m_674139435811179689m_848226107188727543gmail_signature"><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.caam.rice.edu/~mk51/" target="_blank">http://www.caam.rice.edu/~mk51<wbr>/</a><br></div></div></div>
</font></span></font></span></font></span></div></div><span class="gmail-HOEnZb"><font color="#888888">
</font></span></blockquote></div><span class="gmail-HOEnZb"><font color="#888888"><br></font></span></div><span class="gmail-HOEnZb"><font color="#888888">
</font></span></blockquote></div></div></div><span class="gmail-HOEnZb"><font color="#888888"><div><div class="gmail-m_-4173110423624862420m_697185988798551954h5"><br><br clear="all"><div><br></div>-- <br><div class="gmail-m_-4173110423624862420m_697185988798551954m_-4157210347045451781gmail_signature"><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.caam.rice.edu/~mk51/" target="_blank">http://www.caam.rice.edu/~mk51<wbr>/</a><br></div></div></div>
</div></div></font></span></div></div><span class="gmail-HOEnZb"><font color="#888888">
</font></span></blockquote></div><span class="gmail-HOEnZb"><font color="#888888"><br></font></span></div>
</div></div></blockquote></div><br></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature"><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.caam.rice.edu/~mk51/" target="_blank">http://www.caam.rice.edu/~mk51/</a><br></div></div></div>
</div></div>