<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=us-ascii">
<style type="text/css" style="display:none;"><!-- P {margin-top:0;margin-bottom:0;} --></style>
</head>
<body dir="ltr">
<div id="divtagdefaultwrapper" style="font-size:12pt;color:#000000;font-family:Calibri,Helvetica,sans-serif;" dir="ltr">
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Hi Sandeep,</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">intp_rstd(....,1)</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">does not interpolate back.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">It takes the transpose of the forward interpolation operator.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">WHY?</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Because if I want to integrate two functions using quadrature on</p>
<p style="margin-top:0;margin-bottom:0">a fine mesh, I need a routine that interpolates _each_ function to the</p>
<p style="margin-top:0;margin-bottom:0">fine mesh.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Thus, to evaluate</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">     (v,u) := \int  v u dx   =  v^T B u</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">I would use a quadrature rule on M gauss points, where M > N (typically).</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">If u = [ u0 u1 ... uN ]^T is the set of basis coefficients for u, then evaluating</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">u(x) on points xi_k, with xi_k, k=1,...,M being M gauss points, requires:</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">      u(xi_k) = \sum_j=0^N    h_j (xi_k) u_j</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Let uf = [ u(xi_1)  u(xi_2) ... u(xi_M) ]^T be the vector of u values on the</p>
<p style="margin-top:0;margin-bottom:0">fine points.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Then:</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">       uf = J u</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">where J_{kj} :=  h_j(xi_k)</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">and h_j(x) is the jth cardinal Lagrange polynomial of degree N on the</p>
<p style="margin-top:0;margin-bottom:0">spectral element nodal points.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Let  rho_k, k=1,...,M  be the quadrature weights for the M-point gauss rule.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Then,</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">(v,u) ==   \sum_k  vf_k  rho_k  uf_k</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">          =    (Jv)^T  D  Ju</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">          =    v^T  J^T D  J u</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">          =    v^T     B        u</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">where B:=J^T D J, and D = diag (rho_k), k=1,...,M.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">So, when evaluating a matrix-vector product of the form</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">            w = B u</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">we could (and often do) express this as:</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">             w1  =    Ju</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">              w2  =   D w1</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">              w  = J^T  w1</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Now, done.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">intp_rstd(....,1)   implements J^T.</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">hth,</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">Paul</p>
<p style="margin-top:0;margin-bottom:0"><br>
</p>
<p style="margin-top:0;margin-bottom:0">           </p>
</div>
<hr style="display:inline-block;width:98%" tabindex="-1">
<div id="divRplyFwdMsg" dir="ltr"><font face="Calibri, sans-serif" style="font-size:11pt" color="#000000"><b>From:</b> Nek5000-users <nek5000-users-bounces@lists.mcs.anl.gov> on behalf of nek5000-users@lists.mcs.anl.gov <nek5000-users@lists.mcs.anl.gov><br>
<b>Sent:</b> Monday, June 4, 2018 9:24:20 AM<br>
<b>To:</b> nek5000-users@lists.mcs.anl.gov<br>
<b>Subject:</b> [Nek5000-users] intp_rstd; coarse to fine and back</font>
<div> </div>
</div>
<meta content="text/html; charset=utf-8">
<div>
<div dir="ltr">
<div class="x_gmail_default" style="font-size:small">
<div class="x_gmail_default" style="font-size:small">Hello,</div>
<div class="x_gmail_default" style="font-size:small">To know more about the subroutines used for dealiasing, I perform a simple exercise.</div>
<div class="x_gmail_default" style="font-size:small">I take a simple sine function and project it onto the fine mesh and back to the coarse using intp_rstd.</div>
<div class="x_gmail_default" style="font-size:small">When I compare the values of the processed array with the original array, it is scaled by a factor of 1.25204, ie, I get a function like 1.25204*sin(2*x)<br>
</div>
<div class="x_gmail_default" style="font-size:small">Could you please tell me which factor am I missing? Which jacobian has to be multiplied?</div>
<div class="x_gmail_default" style="font-size:small">Attached is the usr file and below are the operations I perform in short.</div>
<div class="x_gmail_default" style="font-size:small"><br>
</div>
<div class="x_gmail_default">      sinax(i,1,1,1)      = sin(2.*x)<br>
</div>
<div class="x_gmail_default">
<div class="x_gmail_default">      </div>
<div class="x_gmail_default">     call intp_rstd(uf1,u1,nx1,nxd,if3d,<wbr>0) ! 0 --> forward</div>
<div class="x_gmail_default">      nxyzd = lxd*lyd*lzd</div>
<div class="x_gmail_default">         do i=1,nxyzd</div>
<div class="x_gmail_default">            tr1   = uf1(i)</div>
<div class="x_gmail_default">           <span> </span><span style="color:rgb(34,34,34); font-family:arial,sans-serif; font-size:small; font-style:normal; font-weight:400; letter-spacing:normal; text-align:start; text-indent:0px; text-transform:none; white-space:normal; word-spacing:0px; background-color:rgb(255,255,255); float:none; display:inline">wf(i)</span><span> </span>=
 tr1*rx(i,1,e)+tr1*rx(i,2,e)+<wbr>tr1*rx(i,3,e)</div>
<div class="x_gmail_default">         enddo<br>
</div>
<div class="x_gmail_default"><br>
</div>
<div class="x_gmail_default">      call intp_rstd(w1,wf,nx1,nxd,if3d,<wbr>1) ! 1 --> back to coarse</div>
<div class="x_gmail_default"><br>
</div>
<div class="x_gmail_default">      nxyz = lx1*ly1*lz1           !              </div>
<div class="x_gmail_default">         do i=1,nxyz               ! to physical space.</div>
<div class="x_gmail_default">            bi = 1./((bm1(i,1,1,e)))   !</div>
<div class="x_gmail_default">            w1(i) = bi*w1(i)</div>
<div class="x_gmail_default">
<div class="x_gmail_default">            jacm=jacmi(i,e) ! I dont know whether to use jacmi(i,e) or not</div>
<div>         enddo</div>
<div><span style="color:rgb(34,34,34); font-family:arial,sans-serif; font-size:12.8px; font-style:normal; font-weight:400; letter-spacing:normal; text-align:start; text-indent:0px; text-transform:none; white-space:normal; word-spacing:0px; background-color:rgb(255,255,255); float:none; display:inline">Thank
 you in advance.</span></div>
<div><span style="color:rgb(34,34,34); font-family:arial,sans-serif; font-size:12.8px; font-style:normal; font-weight:400; letter-spacing:normal; text-align:start; text-indent:0px; text-transform:none; white-space:normal; word-spacing:0px; background-color:rgb(255,255,255); float:none; display:inline"><br>
</span></div>
<div><span style="color:rgb(34,34,34); font-family:arial,sans-serif; font-size:12.8px; font-style:normal; font-weight:400; letter-spacing:normal; text-align:start; text-indent:0px; text-transform:none; white-space:normal; word-spacing:0px; background-color:rgb(255,255,255); float:none; display:inline">Sandeep</span></div>
</div>
</div>
<br clear="all">
</div>
<div>
<div class="x_gmail_signature">
<div dir="ltr">
<div dir="ltr"><br>
</div>
<div>
<div style="margin:2px 0px 0px; font-size:12.8px"></div>
</div>
</div>
</div>
</div>
</div>
</div>
</body>
</html>