<div dir="ltr"><div dir="ltr">On Thu, Mar 11, 2021 at 8:49 AM Mathieu Dutour <<a href="mailto:mathieu.dutour@gmail.com">mathieu.dutour@gmail.com</a>> wrote:<br></div><div class="gmail_quote"><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 dir="ltr">On Thu, 11 Mar 2021 at 13:52, Mark Adams <<a href="mailto:mfadams@lbl.gov" target="_blank">mfadams@lbl.gov</a>> wrote:<br></div><div class="gmail_quote"><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">Mathieu, <div>We have "FieldSplit" support for fields, but I don't know if it has ever been pushed to 1000's of fields so it might fall down. It might work.</div><div>FieldSplit lets you manipulate the ordering, say field major (j) or node major (i).</div></div></blockquote><div>I just looked at it and FieldSplit appears to be used in preconditioner so not exactly relevant.</div><div><br></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>What was unsatisfactory?</div><div>It sounds like you made a rectangular matrix A(1000,3e5) . Is that correct?</div></div></blockquote><div>That is incorrect. The matrix is of size (N, N) with N = 1000 * 3e^5. It is a square</div><div>matrix coming from an implicit scheme.</div><div><br></div><div>Since the other answer appears to have the same misunderstanding, let me try</div><div>to re-explain my point:</div><div>--- In many contexts we need a partial differential equation that is not scalar.</div><div>For example, the shallow water equation has b = 3 fields: H, HU, HV. There are other</div><div>examples like wave modelling where we have something like b = 1000 fields (in a</div><div>discretization).</div><div>--- So, if we work with say an unstructured grid with N nodes then the total number</div><div>of variables of the system will be N_tot = 3N or N_tot = 1000N.</div><div><br></div><div>The linear system has N_tot unknowns and N_tot equations. The entries </div><div>can be written as idx = (i , j) with 1 <= i <= b and 1 <= j <= N.</div></div></div></blockquote><div><br></div><div>This notation does not make sense to me. You are labeling an unknown with an index (i, j) where i in [1, b] and j in [1, N].</div><div>Usually, we label an unknown in a matrix using (i, j) where i, j in [1, N]. Perhaps what you want is a multiindex</div><div><br></div><div> alpha = (i, k) where i in [1, N] and k in [1, b]</div><div><br></div><div>Then a matrix entry is indexed using</div><div><br></div><div> (alpha, beta) = ((i, k), (j, l))</div><div><br></div><div>As Pierre points out, this is exactly what we do in MATBAIJ, where we call "i" a "block index". If you only have block structure,</div><div>then BAIJ is the best you can do. If you have product structure, then KAIJ is better.</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 class="gmail_quote"><div>Thus the non-zero entries in the matrix will be of two kinds:</div><div>--- (idx1, idx2) with idx1 = (i , j) and idx2 = (i' , j) , 1 <= i, i' <= b and 1 <= j <= N.</div><div>Together those define a block in the matrix.</div><div><br></div><div>--- (idx1, idx2) with idx1 = (i , j) and idx2 = (i, j'), 1<= i <= b and 1<= j, j' <= N.</div><div>For each unknown idx1, there will be about 6 unknowns idx2 of this form.</div><div><br></div><div>Otherwise, the block matrices do not have the same coefficients, so a tensor</div><div>product approach does not appear to be workable.</div><div><br></div><div> Mathieu</div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
</blockquote></div>
</blockquote></div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><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.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>