# Inquiries about solving the N-S equation, projection algorithm and parallel

**URL:** https://community.freefem.org/t/inquiries-about-solving-the-n-s-equation-projection-algorithm-and-parallel/3006
**Category:** General Discussion
**Created:** [February 28, 2024, 6:08pm UTC](https://community.freefem.org/t/inquiries-about-solving-the-n-s-equation-projection-algorithm-and-parallel/3006 "2024-02-28T18:08:38Z")
**Posts on this page:** 1
**Showing post:** 9

<div class="post-metadata">

### Author: ![quentin](https://avatars.discourse-cdn.com/v4/letter/q/e47c2d/32.png) [@quentin](https://community.freefem.org/u/quentin)
#### Post date: [April 29, 2024, 6:00am UTC](https://community.freefem.org/t/inquiries-about-solving-the-n-s-equation-projection-algorithm-and-parallel/3006/9 "2024-04-29T06:00:13Z")

</div>

I have written a FreeFEM code to compute vorticity formulation of incompressible flow using `problem`, which is a bit slow and can’t be parallelized.  
So I tried to write a matrix form but met some problem.  
I found an algorithm to solve this problem(JIAN-GUO LIU AND WEINAN E, 2000"SIMPLE FINITE ELEMENT METHOD IN VORTICITY FORMULATION FOR INCOMPRESSIBLE FLOWS"), which decouples the vorticity \omega and stream function \psi :

 ![image](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/5/54e1b4489d994dc7171f25a9b666c6dc70c60ac5.png)  
 ![image](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/c/c99aea639b829af70869078cc6ff85e21389aada.png)  
 ![image](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/6/61892598ad560b6ab0712874d5b9f108542675d2.png)  
 ![image](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/2/25669e6ebb8b626db0256b21aa5a696332dd1d5f.png)  
The X\_{0,h}^k is the subspace of X\_h^k with zero boundary value.  
To compute an auxiliary term \bar{\omega}^{n+1}, I need to get a submatrix (A\_{0,0}A\_{0,b}) from the stiffness matrix A to compute right hand side of (I’), and compute the nonlinear term N.  
In my code using `problem`, the explicit nonlinear term (\nabla \varphi,\omega \nabla^{\perp}\psi) is:

```auto
grad(v1)'*[-dy(psiold),dx(psiold)]*omegaold

```

How to write a `varf` of nonlinear term as (2.6)?  
I know that I can use `varf Stiffness(u,v)=int2d(Th)(grad(u)'*grad(v));` and set `tgv=-10` to build a matrix of Stiffness with zero rows of label 1 , but I have no idea to build A\_{0,b}.  
I found the method to get list of DoF boundary in [Frédéric Hecht](https://community.freefem.org/u/frederichecht)’s reply:

> [@How do I obtain the matrix ordering of the boundaries?](https://community.freefem.org/t/how-do-i-obtain-the-matrix-ordering-of-the-boundaries/2834/2):
>
> I think , you need the list of DoF boundary. border C1(t=0,2\*pi) {x=cos(t); y=sin(t);label=1;} mesh Th=buildmesh(C1(6)); fespace Vh(Th,P1); Vh phi, w; varf vBord(phi,w) = int1d(Th)(w); Vh on1=vBord(0,Vh); for(int i=0; i\< Vh.ndof;++i) if (on1[][i]\>0) cout \<\< i \<\< " " \<\< endl;

If I get the list of DoF on boundary, should I permutate the rows and columns of `Stiffness` to get A？Or is there a better way to implement such algorithm?  
Could you give me some advice?

---

_[View the full topic](https://community.freefem.org/t/inquiries-about-solving-the-n-s-equation-projection-algorithm-and-parallel/3006)._