Global computation of an error with PETSc using block matrices

Hello,
I’ve adapted the part of the stokes-block-2d-PETSc,edp to conduct a convergence test with a manufactured solution.
I need to have tgv=-1. The code works out OK, including the plots, but the errors are almost stagnant.
I would appreciate it if anyone could point out the issue in the code. Thanks in advance.

MWE.edp (4.9 KB)

--- wrong.edp	2023-06-19 06:30:30
+++ right.edp	2023-06-19 06:30:26
@@ -90,8 +90,9 @@
 dC = vPb22(VhP, VhP);
 matrix<real> B = vPb12(VhP, VhV);
 Mat dB(dA, dC, B);
+Mat dD(dA, dC, B);
 {
-  varf vonBorder([u, uB], [v, vB]) = on(1, u = u1exact, uB = u2exact);
+  varf vonBorder([u, uB], [v, vB]) = on(1, u = 1, uB = 1);
   real[int] onBorder = vonBorder(0, VhV, tgv=-1);
   real[int] onBorderPETSc;
   ChangeNumbering(dA, onBorder, onBorderPETSc);
@@ -100,7 +101,7 @@
 
 Mat dBT(dA, dC, B);
 Mat dS = [[dA , dB],
-	  [dBT', dC ]];
+	  [dD', dC ]];
 real[int] rhsPETSc;
 {
   real[int] rhsV = vPb11(0, VhV, tgv = -1);
@@ -152,5 +153,5 @@
 if(mpirank == 0){
   savevtk("Error.vtu", ThGlobal, [solVGlobal,solVGlobalB],solPGlobal,communicator = mpiCommSelf);
   cout.scientific << "L2 v error = " << sqrt(int2d(ThGlobal)(solVGlobal^2+solVGlobalB^2)/ThGlobal.measure) << endl;
-  cout.scientific << "L2 p error = " << sqrt(int2d(ThGlobal)(sumP^2)/ThGlobal.measure) << endl;
+  cout.scientific << "L2 p error = " << sqrt(int2d(ThGlobal)(solPGlobal^2)/ThGlobal.measure) << endl;
  }
1 Like

Thank you!
I see the problem now.