# Poor results with PETSc

**URL:** <https://community.freefem.org/t/poor-results-with-petsc/1058>\
**Category:** General Discussion\
**Created:** [June 24, 2021, 4:24pm UTC](https://community.freefem.org/t/poor-results-with-petsc/1058 "2021-06-24T16:24:29Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![mazzman](https://avatars.discourse-cdn.com/v4/letter/m/b9bd4f/32.png) [@mazzman](https://community.freefem.org/u/mazzman)\
**Post date:** [June 24, 2021, 4:24pm UTC](https://community.freefem.org/t/poor-results-with-petsc/1058/1 "2021-06-24T16:24:29Z")

</div>

Hello,  
I implemented a transient 3D piezoanalysis that works fine with `problem`. Now, I want to accelerate the simulation by using a parallel solver. However, I can’t get proper results using PETSc. The results seem to not converge properly, as I keep getting a very distorted mesh. I tried different solvers like lu or gmres. Also a “static” simulation (finding displacement U) does converge okay’isch.

Only the results in the first time step are accuarte (with lu). Could someone please help me with this problem?

```auto
load "msh3"
load "lapack"
load "PETSc"

//...

varf piezoTrans([phi, Fx,Fy,Fz], [Vphi, Vx, Vy, Vz]) =
	 int3d(Th)( dings(sigma(Ux,Uy,Uz), [Vx,Vy,Vz])+ dings2(D1(Ux,Uy,Uz),Vphi) )  
	 +int3d(Th)( dings2(D2(phi),Vphi) - dings(sigmaE(phi), [Vx,Vy,Vz]))
			  
	
	+int3d(Th)(Fx*Vx + Fy*Vy + Fz*Vz)
	
	 +on(elektrode1, phi=50)
	 +on(elektrode2, phi=0);

//this is exactly the same as above, but works
Vh phi2, Fx2,Fy2,Fz2, Vphi, Vx, Vy, Vz;
problem piezoTrans2([phi2, Fx2,Fy2,Fz2], [Vphi, Vx, Vy, Vz]) =
	 int3d(Th)( dings(sigma(Ux,Uy,Uz), [Vx,Vy,Vz])+ dings2(D1(Ux,Uy,Uz),Vphi) )  
	 +int3d(Th)( dings2(D2(phi2),Vphi) - dings(sigmaE(phi2), [Vx,Vy,Vz]))
			  
	
	+int3d(Th)(Fx2*Vx + Fy2*Vy + Fz2*Vz)
	
	 +on(elektrode1, phi2=50)
	 +on(elektrode2, phi2=0);

fespace Vh2(Th, [P1,P1,P1,P1]);

real dt = 20e-9;
real T = 10e-6;
int m, M = T/dt;

Vh Ux, Uy, Uz;

mesh3[int] anim(M);

matrix AA = piezoTrans(Vh2,Vh2, tgv=-1);
real[int] rhs = piezoTrans(0,Vh2, tgv=-1);
Mat A(AA);

for(m = 0; m < 3; m++){	
	AA = piezoTrans(Vh2,Vh2, tgv=-1);
	A = Mat(A);
	rhs = piezoTrans(0,Vh2, tgv=-1);
	
	
	
	set(A, sparams="-pc_type lu -ksp_view -ksp_monitor -ksp_monitor_true_residual -ksp_view_final_residual");
	//set(A, sparams="-ksp_view -ksp_monitor -ksp_monitor_true_residual -ksp_view_final_residual -pc_type gamg -ksp_type gmres");
	//set(A, sparams="-ksp_view -ksp_monitor -ksp_monitor_true_residual -ksp_view_final_residual -pc_type lu -ksp_type preonly -ksp_max_it 200 -n 4");
	
	Vh2 [phi,Fx,Fy,Fz];
	phi[] = A^-1*rhs;
	
	piezoTrans2;
	
	Vh deltaX, deltaY, deltaZ;
	
	deltaX = Fx-Fx2;
	deltaY = Fy-Fy2;
	deltaZ = Fz-Fz2;
	
	cout << "AbsX: " << (deltaX[]'*deltaX[])/deltaX[].n << endl;
	cout << "AbsY: " << (deltaY[]'*deltaY[])/deltaY[].n << endl;
	cout << "AbsZ: " << (deltaZ[]'*deltaZ[])/deltaZ[].n << endl;
	cout << "Rel: " << (Fx(r2,0,0)-Fx2(r2,0,0))/Fx2(r2,0,0) << endl;
	
	//v1x = Fx/rho*dt;
	//v1y = Fy/rho*dt;
	//v1z = Fz/rho*dt;
	
	v1x = v1x + Fx2/rho*dt;
	v1y = v1y + Fy2/rho*dt;
	v1z = v1z + Fz2/rho*dt;
	
	Ux = Ux + v1x * dt;
	Uy = Uy + v1y * dt;
	Uz = Uz + v1z * dt;
	

	//just some animation stuff here
	//...
} 

```

---

<div class="post-metadata">

**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\
**Post date:** [June 24, 2021, 4:42pm UTC](https://community.freefem.org/t/poor-results-with-petsc/1058/2 "2021-06-24T16:42:51Z")

</div>

You have to be careful when switching from `problem` to `varf`, you have to compute the opposite of your linear form. Have you checked that? `piezoTrans` and `piezoTrans2` seems to be exactly the same, so any code using `varf` (you can check without PETSc) will not work.

---

<div class="post-metadata">

**Author:** ![mazzman](https://avatars.discourse-cdn.com/v4/letter/m/b9bd4f/32.png) [@mazzman](https://community.freefem.org/u/mazzman)\
**Post date:** [June 24, 2021, 5:02pm UTC](https://community.freefem.org/t/poor-results-with-petsc/1058/3 "2021-06-24T17:02:24Z")

</div>

Sorry, I don’t really understand what you mean by “the opposite of my linear form”.

---

<div class="post-metadata">

**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\
**Post date:** [June 24, 2021, 5:19pm UTC](https://community.freefem.org/t/poor-results-with-petsc/1058/4 "2021-06-24T17:19:26Z")

</div>

Multiply by -1. See [FreeFem-sources/Laplace.edp at develop · FreeFem/FreeFem-sources · GitHub](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/tutorial/Laplace.edp#L9) vs. [FreeFem-sources/first.edp at develop · FreeFem/FreeFem-sources · GitHub](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/3d/first.edp#L93).

---

<div class="post-metadata">

**Author:** ![mazzman](https://avatars.discourse-cdn.com/v4/letter/m/b9bd4f/32.png) [@mazzman](https://community.freefem.org/u/mazzman)\
**Post date:** [June 24, 2021, 5:53pm UTC](https://community.freefem.org/t/poor-results-with-petsc/1058/5 "2021-06-24T17:53:50Z")

</div>

Thank you very much! That finally solved my problem!
