<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Fri, Jun 1, 2018 at 12:07 PM, Samuel Lanthaler <span dir="ltr"><<a href="mailto:s.lanthaler@gmail.com" target="_blank">s.lanthaler@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
  
    
  
  <div bgcolor="#FFFFFF" text="#000000">
    <div class="m_3362958314088522262moz-cite-prefix">On 06/01/2018 03:42 PM, Matthew Knepley
      wrote:<br>
    </div>
    <blockquote type="cite">
      <div dir="ltr">
        <div class="gmail_extra">
          <div class="gmail_quote">On Fri, Jun 1, 2018 at 9:21 AM,
            Samuel Lanthaler <span dir="ltr"><<a href="mailto:s.lanthaler@gmail.com" target="_blank">s.lanthaler@gmail.com</a>></span>
            wrote:<br>
            <blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hi,<br>
              <br>
              I was wondering what the most efficient way to use MatPtAP
              would be in the following situation: I am discretizing a
              PDE system. The discretization yields a matrix A that has
              a band structure (with k upper and lower bands, say). In
              order to implement the boundary conditions, I use a
              transformation matrix P which is essentially the unit
              matrix, except for the entries P_{ij} where i,j<k and
              n-i,n-j<k, so<br>
              <br>
              P =  [ B, 0, 0, 0, ..., 0, 0 ]<br>
                     [  0, 1, 0, 0, ..., 0, 0 ]<br>
                     [                              ]<br>
                     [                              ]<br>
                     [                  ..., 1, 0 ]<br>
                     [  0, 0, 0, 0, ..., 0, C ]<br>
              <br>
              with B,C are (k-by-k) matrices.<br>
              Right now, I'm simply constructing A, P and calling<br>
              <br>
                  CALL MatPtAP(petsc_matA,petsc_matP,<wbr>MAT_INITIAL_MATRIX,PETSC_DEFAU<wbr>LT_REAL,petsc_matPtAP,ierr)<br>
              <br>
              where I haven't done anything to pestc_matPtAP, prior to
              this call. Is this the way to do it?<br>
              <br>
              I'm asking because, currently, setting up the matrices A
              and P takes very little time, whereas the operation
              MatPtAP is taking quite long, which seems very odd... The
              matrices are of type MPIAIJ. In my problem, the total
              matrix dimension is around 10'000 and the matrix blocks
              (B,C) are of size ~100.<br>
            </blockquote>
            <div><br>
            </div>
            <div>Are you sure this is what you want to do? Usually BC
              are local, since by definition PDE are local, and</div>
            <div>are applied pointwise. What kind of BC do you have
              here?</div>
            <div><br>
            </div>
          </div>
        </div>
      </div>
    </blockquote>
    <br>
    The boundary conditions are a mixture of Dirichlet and Neumann; in
    my case, the PDE is a system involving 8 variables on a disk, where
    the periodic direction is discretized using a Fourier series
    expansion, the radial direction uses B-splines. <br>
    <br>
    In reality, I have two matrices A,B, and want to solve the
    eigenvalue problem \lambda*B*x = A*x. <br>
    I found it quite convenient to use a transformation P to a different
    set of variables y, such that x=P*y and x satisfies the BC iff
    certain components of y are 0. The latter is enforced by inserting
    spurious eigenvalues at the relevant components of y in the
    transformed eigenvalue problem \lambda*Pt*B*P*y=Pt*A*P*y. After
    solving the EVP in terms of y, I get back x=P*y.<br>
    Is this an inherently bad/inefficient way of enforcing BC's? Thanks.<br></div></blockquote><div><br></div><div>I cannot quite see what is happening. However it appears that you are coupling all the unknowns on a boundary to enforce the boundary condition.</div><div>If this is true, its pretty expensive. For example, in 2D you have N^2 unknowns for N boundary unknowns. However, dense coupling on the boundary</div><div>produces N^2 storage and N^3 computation, which is not great. Am I missing something? Its fine to do this for some sizes (2D BEM guys do it all</div><div>the time), but for bigger stuff it gets untenable.</div><div><br></div><div>  Thanks,</div><div><br></div><div>    Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div bgcolor="#FFFFFF" text="#000000">
    <br>
    <br>
    <blockquote type="cite">
      <div dir="ltr">
        <div class="gmail_extra">
          <div class="gmail_quote">
            <div>  Thanks,</div>
            <div><br>
            </div>
            <div>    Matt</div>
            <div> </div>
            <blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
              Thanks in advance for any ideas.<br>
              <br>
              Cheers,<br>
              Samuel<br>
            </blockquote>
          </div>
          <br>
          <br clear="all"><span class="HOEnZb"><font color="#888888">
          <div><br>
          </div>
          -- <br>
          <div class="m_3362958314088522262gmail_signature" data-smartmail="gmail_signature">
            <div dir="ltr">
              <div>
                <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/%7Emk51/" target="_blank">https://www.cse.buffalo.edu/~<wbr>knepley/</a><br>
                  </div>
                </div>
              </div>
            </div>
          </div>
        </font></span></div>
      </div>
    </blockquote>
    <br>
  </div>

</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div><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">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div>
</div></div>