# Wrong answer from imposing a point value for the pressure via the penalty method

**URL:** https://community.freefem.org/t/wrong-answer-from-imposing-a-point-value-for-the-pressure-via-the-penalty-method/3176
**Category:** General Discussion
**Created:** [April 30, 2024, 2:58am UTC](https://community.freefem.org/t/wrong-answer-from-imposing-a-point-value-for-the-pressure-via-the-penalty-method/3176 "2024-04-30T02:58:15Z")
**Posts on this page:** 1
**Page:** 1

<div class="post-metadata">

### Author: ![ars\_cn](https://avatars.discourse-cdn.com/v4/letter/a/5f9b8f/32.png) [@ars\_cn](https://community.freefem.org/u/ars_cn)
#### Post date: [April 30, 2024, 2:58am UTC](https://community.freefem.org/t/wrong-answer-from-imposing-a-point-value-for-the-pressure-via-the-penalty-method/3176/1 "2024-04-30T02:58:15Z")

</div>

Hello,

I am working on the Stokes equation now. I noticed FreeFEM introduces an extra small L2 term on the pressure to stabilize the scheme. However, I am trying to impose the pressure with a point value to gain the well-posedness of the system. I found the d.o.f of pressure via the following approach (sorry, the forum website does not allow me to upload the file, so I will paste my code here):

//\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*  
fespace Xh(Th,[P2,P2,P1]);  
Xh [uind,vind,pind] = [1,2,3];  
for(int irow=0;irow\<Ahat.n;irow++){[//Xh.ndof](https://Xh.ndof)  
//int irow = 2;  
if(uind[irow]==3){  
/\* Change to p = pvalue!  
irowInd = irow;  
//cout\<\<"Ahat = "\<\<Ahat(irow,irow)\<\<endl;  
for(int icol=0;icol\<Ahat.n;icol++){  
if(abs(Ahat(irow,icol))\>1e-8){  
Ahat(irow,icol) = 0.0;  
}  
}  
for(int icol=0;icol\<Ahat.n;icol++){  
if(abs(Ahat(icol,irow))\>1e-8){  
bhat(icol) = bhat(icol) - Ahat(icol,irow)_ppiont(nu, t,ux[][irow],uy[][irow]);  
}  
}  
Ahat(irow,irow) = 1.0;  
bhat(irow) = ppiont(nu, t,ux[][irow],uy[][irow]);//tgv_  
_/  
// Penalty approach  
Ahat(irow,irow) = 1e30;  
bhat(irow) = 1e30_ppiont(nu, t,ux[irow],uy[irow]);  
break;  
}  
}[//Xh.ndof](https://Xh.ndof)  
solhat = Ahat^-1_bhat;  
uhat[] = solhat;  
//_\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*\*  
Unfortunately, it does not give the exact solution (ppiont(nu, t,ux[irow],uy[irow])) at uhat[irowInd]. As you can see, I also try the other way to impose such boundary condition, but it also fails. Can anyone help me out?
