@@ -167,40 +167,42 @@ int main(int argc, char *argv[])
167167 gt1->UseEA (use_ea);
168168 gt2->UseEA (use_ea);
169169
170- const Operator &R1 = gt1->ForwardOperator (); // Restriction
171- const Operator &R2 = gt2->ForwardOperator (); // Restriction
170+ const Operator &P1 = gt1->BackwardOperator (); // Prolongation 1 (LO->IM)
171+ const Operator &P2 = gt2->ForwardOperator (); // Prolongation 2 (IM->HO)
172172
173- const Operator &P1 = gt1->BackwardOperator (); // Prolongation
174- const Operator &P2 = gt2->BackwardOperator (); // Prolongation
173+ const Operator &R1 = gt1->ForwardOperator (); // Restriction 1 (IM->LO)
174+ const Operator &R2 = gt2->BackwardOperator (); // Restriction 2 (HO->IM)
175175
176176 // ======================================================
177177 // Compute GridFunctions
178178 // ======================================================
179179
180180 std::cout<<" SK: LO projection" <<std::endl;
181181
182- // STEP1:LO projection
183- FunctionCoefficient RHO (RHO_exact);
182+ // STEP1: LO projection
184183
184+ FunctionCoefficient RHO (RHO_exact);
185185 rho_lo.ProjectCoefficient (RHO);
186186 rho_lo.SetTrueVector ();
187187 rho_lo.SetFromTrueVector ();
188188
189189 std::cout<<" SK: LO compute mass" <<std::endl;
190190
191- real_t mass_lo = compute_mass (&fespace_lo, -1.0 , dc_lo, " LO " );
191+ real_t mass_lo = compute_mass (&fespace_lo, -1.0 , dc_lo, " LO" );
192192 if (vis) { visualize (dc_lo, " LO" , Wx, Wy, visport); Wx += offx; }
193193
194+ // STEP2: Prolongation 1 (LO->IM)
194195
195- // STEP2: LO->IM Prolongation
196- direction = " LO -> IM @ IM" ;
197196 P1.Mult (rho_lo, rho_im); // rho_im = P1 * rho
198197
199- // P1.Mult(rho_lo, rho_im); // rho_im = P1 * rho_lo
200- // real_t mass_im = compute_mass(&fespace_im, -1.0, dc_im, "P(R(HO)) ");
201-
202- // compute_mass(&fespace_lo, mass_im, dc_lo, "R(HO) ");
203- // if (vis) { visualize(dc_lo, "R(HO)", Wx, Wy, visport); Wx += offx; }
198+ real_t mass_im = compute_mass (&fespace_im, -1.0 , dc_im, " P1(LO) " );
199+ if (vis) { visualize (dc_im, " IM=P1(LO)" , Wx, Wy, visport); Wx += offx; }
200+
201+ // STEP3: Prolongation 2 (IM->HO)
202+
203+ P2.Mult (rho_im, rho_hi); // rho_hi = P2 * rho_im
204+ real_t mass_hi = compute_mass (&fespace_hi, -1.0 , dc_hi, " P2(IM) " );
205+ if (vis) { visualize (dc_hi, " HO=P2(IM)" , Wx, Wy, visport); Wx += offx; }
204206
205207 // if (gt1->SupportsBackwardsOperator()) // 1
206208 // {
0 commit comments