<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=Windows-1252">
</head>
<body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space;">
Hi Paul and Stefan,
<div><br>
</div>
<div>Thank you for the explanations, things are going gradually more clear in my mind.</div>
<div><br>
</div>
<div>So for a practical application, still in the context of cvode, I would like to add a new term to the RHS of the temperature equation. This new term is roughly</div>
<div><br>
</div>
<div>(vdiff * grad (Y)) grad (T) </div>
<div>where Y is a passive scalar, vdiff the diffusion coefficient of this passive scalar and T the temperature.</div>
<div><br>
</div>
<div>Basically in the routine makeq.f the implementation should looks like to</div>
<div><br>
</div>
<div>call wgradm1(Y_x, Y_y, Y_z, t(1,1,1,1,ifield+1), nelv)</div>
<div>call dssum(grad_Y, Y_x, Y_y, Y_z)</div>
<div>call col2(grad_Y,bintm1,ntot)</div>
<div><br>
</div>
<div>call col2(grad_Y, vdiff(1,1,1,1,ifield+2) ,ntot )</div>
<div><br>
</div>
<div>
<div>call wgradm1(T_x, T_y, T_z, t(1,1,1,1,ifield), nelv)</div>
<div>call dssum(grad_T, T_x, T_y, T_z)</div>
</div>
<div><br>
</div>
<div>call col2(grad_T,grad_Y,ntot)</div>
<div><br>
</div>
<div>call add2(bq(1,1,1,1,ifield), grad_T,ntot)</div>
<div><br>
</div>
<div>and then convab and wlaplacian routines will add the other terms, and it will finish with </div>
<div>
<div>call dssum(bq,nx1,ny1,nz1)</div>
<div> call col2(bq,bintm1,ntot)</div>
<div> call col2(bq,bm1,ntot)</div>
</div>
<div><br>
</div>
<div>Am I right??</div>
<div>Thank you for your help !</div>
<div><br>
</div>
<div><br>
<div>
<div style="color: rgb(0, 0, 0); letter-spacing: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px; word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space;">
<div>Best regards,</div>
<div>Emmanuel</div>
</div>
<br class="Apple-interchange-newline">
<br class="Apple-interchange-newline">
</div>
<br>
<div>
<div>On 5 Jun 2014, at 19:25, <<a href="mailto:nek5000-users@lists.mcs.anl.gov">nek5000-users@lists.mcs.anl.gov</a>> <<a href="mailto:nek5000-users@lists.mcs.anl.gov">nek5000-users@lists.mcs.anl.gov</a>> wrote:</div>
<br class="Apple-interchange-newline">
<blockquote type="cite">
<div style="font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px; direction: ltr; font-family: Tahoma; font-size: 10pt;">
<br class="Apple-interchange-newline">
Hi Emmanuel,
<div><br>
</div>
<div>Note that bintm1 and bm1 are not precise inverses of one another.</div>
<div>bm1 is the unassembled mass matrix - i.e., it is the collection of local</div>
<div>(to each element) mass matrices, each of which is diagonal.  </div>
<div><br>
</div>
<div>
<div>Conceptually, we denote the local mass matrix<span style="font-family: Tahoma; font-size: 10pt;"> as B_L, for local.  </span></div>
<div><span style="font-family: Tahoma; font-size: 10pt;">Then, the assembled mass matrix is:</span></div>
<div><br>
</div>
<div>     B := Q^T B_L Q</div>
<div><br>
</div>
<div>where Q is the Boolean matrix that copies global degrees of freedom</div>
<div>to their shared counterparts (i.e., copies a values to shared faces, </div>
<div>edges, vertices, etc.).   For the SEM, B_L is diagonal, as is B, so</div>
<div>Binv := B^-1 is easily computed.   Once the entries of this diagonal</div>
<div>matrix are known, we can copy them to their element-based locations,</div>
<div>so that we can generate B^-1 * r very easily.  </div>
<div><br>
</div>
</div>
<div>Often, products with these matrices, in conjunction with dssum, are used</div>
<div>to generate continuous fields and that appears to be what is happening in </div>
<div>this case.   An example would be if you want to find a the continuous </div>
<div>counterpart to a discontinuous function g.  The usual projection approach</div>
<div>is:  Find gc \in H1 such that  \int v gc dV = \int v g dV for all v in H1, which</div>
<div>ends up looking like:</div>
<div><br>
</div>
<div>Find gc \in H1 such that:</div>
<div><br>
</div>
<div><br>
</div>
<div>         (v,gc) = (v,g)  for all g \in H1</div>
<div><br>
</div>
<div>or</div>
<div><br>
</div>
<div>         v B gc = v Q^T B_L g_L</div>
<div><br>
</div>
<div>               gc = B^-1 Q^T B_L g_L</div>
<div><br>
</div>
<div>               gc_L = Q B^-1 Q^T B_L g_L</div>
<div><br>
</div>
<div>In Nek, the last statement is expressed as:</div>
<div><br>
</div>
<div>              call col2(g,bm1,n)</div>
<div>              call dssum(g,nx1,ny1,nz1)</div>
<div>              call col2(g,bintm1,n)</div>
<div><br>
</div>
<div>In your cvode case, these statements are permuted because of the</div>
<div>initial state of the data and the desired output (which is to be B*q --</div>
<div>which is why the variable is labeled bq).</div>
<div><br>
</div>
<div>Paul</div>
<div><br>
</div>
<div><br>
</div>
<div><br>
</div>
<div><br>
<div style="font-family: 'Times New Roman'; font-size: 16px;">
<hr tabindex="-1">
<div id="divRpF599061" style="direction: ltr;"><font face="Tahoma" size="2"><b>From:</b><span class="Apple-converted-space"> </span><a href="mailto:nek5000-users-bounces@lists.mcs.anl.gov">nek5000-users-bounces@lists.mcs.anl.gov</a><span class="Apple-converted-space"> </span>[<a href="mailto:nek5000-users-bounces@lists.mcs.anl.gov">nek5000-users-bounces@lists.mcs.anl.gov</a>]
 on behalf of<span class="Apple-converted-space"> </span><a href="mailto:nek5000-users@lists.mcs.anl.gov">nek5000-users@lists.mcs.anl.gov</a><span class="Apple-converted-space"> </span>[<a href="mailto:nek5000-users@lists.mcs.anl.gov">nek5000-users@lists.mcs.anl.gov</a>]<br>
<b>Sent:</b><span class="Apple-converted-space"> </span>Thursday, June 05, 2014 2:38 AM<br>
<b>To:</b><span class="Apple-converted-space"> </span><a href="mailto:nek5000-users@lists.mcs.anl.gov">nek5000-users@lists.mcs.anl.gov</a><br>
<b>Subject:</b><span class="Apple-converted-space"> </span>Re: [Nek5000-users] weak derivative, dssum and binvm1 in cvode_driver<br>
</font><br>
</div>
<div></div>
<div>
<div>Hi Emmanuel</div>
<div><br>
</div>
<div>From what I recall: </div>
<div>CVODE solves the governing non-linear eqns in the strong form. That’s why we have to project every term of the RHS onto H2 (dssum and collocate with bintm1). The sum is projected again onto H2 because it gets devided by vtrans. Hope this helps.</div>
<div><br>
</div>
<div>Cheers,</div>
<div>Stefan</div>
<div><br>
</div>
<span id="OLK_SRC_BODY_SECTION">
<div style="font-family: Calibri; font-size: 11pt; text-align: left; border-width: 1pt medium medium; border-style: solid none none; padding: 3pt 0in 0in; border-top-color: rgb(181, 196, 223);">
<span style="font-weight: bold;">Von:<span class="Apple-converted-space"> </span></span><<a href="mailto:nek5000-users@lists.mcs.anl.gov" target="_blank">nek5000-users@lists.mcs.anl.gov</a>><br>
<span style="font-weight: bold;">Antworten an:<span class="Apple-converted-space"> </span></span><<a href="mailto:nek5000-users@lists.mcs.anl.gov" target="_blank">nek5000-users@lists.mcs.anl.gov</a>><br>
<span style="font-weight: bold;">Datum:<span class="Apple-converted-space"> </span></span>Thursday, June 5, 2014 at 9:20 AM<br>
<span style="font-weight: bold;">An:<span class="Apple-converted-space"> </span></span>"<a href="mailto:nek5000-users@lists.mcs.anl.gov" target="_blank">nek5000-users@lists.mcs.anl.gov</a>" <<a href="mailto:nek5000-users@lists.mcs.anl.gov" target="_blank">nek5000-users@lists.mcs.anl.gov</a>><br>
<span style="font-weight: bold;">Betreff:<span class="Apple-converted-space"> </span></span>Re: [Nek5000-users] weak derivative, dssum and binvm1 in cvode_driver<br>
</div>
<div><br>
</div>
<div>
<div style="word-wrap: break-word;">I’m adding another quick question:
<div><br>
</div>
<div>why in the routine makeq, under the cvode statement, there are these following 2 lines:</div>
<div>
<div>             call col2(bq,bintm1,ntot)</div>
<div>             call col2(bq,bm1,ntot)</div>
<div><br>
</div>
<div>I’m not sure to understand why we are multiplying bq by the inverse of the mass matrix (bintm1) and just after multiplying it back by the mass matrix (bm1). Moreover, after leaving the makeq routine, fctfun will multiplying again ydot=bq/vtrans by the
 inverse of the mass matrix.</div>
<div><br>
</div>
<div>Why all of these operations? Thank you for your help !</div>
<div><br>
</div>
<div><br>
</div>
<div><br>
</div>
<div>
<div style="letter-spacing: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; word-wrap: break-word;">
<div>Best regards,</div>
<div>Emmanuel</div>
</div>
<br class="Apple-interchange-newline">
<br class="Apple-interchange-newline">
</div>
<br>
<div>
<div>On 5 Jun 2014, at 12:14, <<a href="mailto:nek5000-users@lists.mcs.anl.gov" target="_blank">nek5000-users@lists.mcs.anl.gov</a>> <<a href="mailto:nek5000-users@lists.mcs.anl.gov" target="_blank">nek5000-users@lists.mcs.anl.gov</a>> wrote:</div>
<br class="Apple-interchange-newline">
<blockquote type="cite">
<div style="word-wrap: break-word;">Dear Nek team,
<div><br>
</div>
<div>I have a quick question: in cvode_driver.f, in the routine fcvfun, do all the derivatives must be in weak formulation (using wgradm1) when constructing the right hand side term ydot?</div>
<div>Because at the end of the routine fcvfun,  dssum and multiplication by binvm1 are applied to the right hand side term ydot.</div>
<div><br>
</div>
<div>It’s also not clear to me if I need to apply dssum and binvm1 for each derivative term appearing in my equation (in other terms, each time after calling wgradm1), or if I need to apply dssum and binvm1 only to the whole RHS term ydot ?</div>
<div><br>
</div>
<div>I’m afraid of inconsistency if apply twice dssum and binvm1: one during the construction of operators and one on the final ydot term.</div>
<div><br>
</div>
<div>Thank you for you help !</div>
<div><br>
<div>
<div style="letter-spacing: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; word-wrap: break-word;">
<div>Best regards,</div>
<div>Emmanuel</div>
</div>
<br class="Apple-interchange-newline">
<br class="Apple-interchange-newline">
</div>
<br>
</div>
</div>
_______________________________________________<br>
Nek5000-users mailing list<br>
<a href="mailto:Nek5000-users@lists.mcs.anl.gov" target="_blank">Nek5000-users@lists.mcs.anl.gov</a><br>
<a href="https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users" target="_blank">https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users</a><br>
</blockquote>
</div>
<br>
</div>
</div>
</div>
_______________________________________________ Nek5000-users mailing list<a href="mailto:Nek5000-users@lists.mcs.anl.gov" target="_blank">Nek5000-users@lists.mcs.anl.gov</a><span class="Apple-converted-space"> </span><a href="https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users" target="_blank">https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users</a></span></div>
</div>
</div>
</div>
<span style="font-family: Arial, sans-serif; font-size: 14px; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px; float: none; display: inline !important;">_______________________________________________</span><br style="font-family: Arial, sans-serif; font-size: 14px; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px;">
<span style="font-family: Arial, sans-serif; font-size: 14px; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px; float: none; display: inline !important;">Nek5000-users
 mailing list</span><br style="font-family: Arial, sans-serif; font-size: 14px; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px;">
<a href="mailto:Nek5000-users@lists.mcs.anl.gov" style="font-family: Arial, sans-serif; font-size: 14px; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px;">Nek5000-users@lists.mcs.anl.gov</a><br style="font-family: Arial, sans-serif; font-size: 14px; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px;">
<a href="https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users" style="font-family: Arial, sans-serif; font-size: 14px; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px;">https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users</a></blockquote>
</div>
<br>
</div>
</body>
</html>