<!DOCTYPE html>
<html>
<head>
<meta http-equiv="content-type" content="text/html; charset=UTF-8">
</head>
<body>
<p>Hi there,</p>
<p>I am (successfully) using the following piece of code to insert
elemental matrices into the "global" matrix "jac":</p>
<p>
<blockquote type="cite">
<span style="font-family:monospace"><span
style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">if(
matrix_assembly )then</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">i = 0</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">PetscCall(DMPlexGetTransitiveClosure(dm,
cell, useCone, PETSC_NULL_INTEGER, p_points, ierr))</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!
loop over the vertices of the cell</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">do j =
vtxbgn, vtxend, 2 ! this is ndim-dependent</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
i = i + 1</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
point = p_points(j)</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">PetscCall(PetscSectionGetDof(section,
point, ndofs, ierr))</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
PetscCall(PetscSectionGetOffset(section,
point, offset, ierr))</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
row(i) = offset/ndofs</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
col(i) = offset/ndofs</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
end do</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">do i = 1,
nofvert</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">do j = 1,
nofvert</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">dummy =
transpose(eltmat(:, :, i, j))</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">ir = row(i)</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">ic = col(j)</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">p_dummy</span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">is</span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">a</span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">1d</span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">view</span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">of</span><span
style="color:#000000;background-color:#ffffff;"> </span><span
style="color:#000000;background-color:#ffffff;">dummy</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">PetscCall(MatSetValuesBlockedLocal(jac,
2, [ir], 1, [ic], p_dummy, ADD_VALUES, ierr))</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">end do</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">end do</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">PetscCall(DMPlexRestoreTransitiveClosure(dm,
cell, useCone, PETSC_NULL_INTEGER, p_points, ierr))</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">!</span><span
style="color:#000000;background-color:#ffffff;">
</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><span
style="color:#000000;background-color:#ffffff;">end if !
matrix_assembly</span><br>
<span style="color:#000000;background-color:#ffffff;">
</span><br>
</span>
</blockquote>
I would however be happy to replace all of the above with a single
call to:</p>
<p>
<blockquote type="cite">
<span style="font-family:monospace"><span
style="color:#000000;background-color:#ffffff;">PetscCall(DMPlexMatSetClosure(dm,
PETSC_NULL_SECTION, PETSC_NULL_SECTION, jac, cell, values,
ADD_VALUES, ierr ))</span><br>
</span>
</blockquote>
which I actually succeeded to do as long as ndofs = 1.</p>
<p>The question is:</p>
<p>how (in which order) should I "transfer" the entries of
eltmat(1:ndofs,1:ndofs,1:nofvert,1:nofvert) into the values array
required by DMPlexMatSetClosure ?</p>
<p>Thanks,</p>
<p>Aldo</p>
<pre class="moz-signature" cols="72">--
Dr. Aldo Bonfiglioli
Associate professor of Fluid Mechanics
Dipartimento di Ingegneria
Universita' della Basilicata
V.le dell'Ateneo Lucano, 10 85100 Potenza ITALY
tel:+39.0971.205203 fax:+39.0971.205215
web: <a class="moz-txt-link-freetext" href="https://urldefense.us/v3/__http://docenti.unibas.it/site/home/docente.html?m=002423__;!!G_uCfscf7eWS!ZusjOf_nCidORTOchtNsX_xGJZNjiAXMKKmn4rDhnR-3Y6mFzMF9Z7f7iJFeGj6tgcBzT2L5gtfxO56AMvPr3pgXHaCPziR-QL4$">http://docenti.unibas.it/site/home/docente.html?m=002423</a></pre>
</body>
</html>