\documentclass{template-jurnal} %Bagian ini diedit oleh editor Jurnal Matematika UNAND \hyphenation{di-tulis-kan di-terang-kan} %Perhatikan aturan penulisan dan ukuran huruf yang digunakan \begin{document} \markboth{Bayu Prihandono et al.} %Jika lebih dari dua penulis, tuliskan sebagai Nama Penulis Pertama dkk. {Existence and uniqueness in the linearised one and two-dimensional problem of partial differential...} %%%%%%%%%%%%%%%%%%%%% Publisher's Area please ignore %%%%%%%%%%%%%%% % \catchline{}{}{}{}{} % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \title{\uppercase{Existence and uniqueness in the linearised one and two-dimensional problem of partial differential equations with variational method}} \author{BAYU PRIHANDONO$^{a}$\footnote{penulis korespondensi}, MARIATUL KIFTIAH$^{b}$, YUDHI$^{c}$\\ (ditulis dengan huruf besar, \textbf{TANPA GELAR})} \address{$^{a}$ Universitas Tanjungpura,\\ $^{b}$ Universitas Tanjungpura,\\ $^{c}$ Universitas Tanjungpura.\\ email : \email{bayuprihandono@math.untan.ac.id, kiftiahmariatul@math.untan.ac.id, yudhi@math.untan.ac.id}} \maketitle \setcounter{page}{1} %bagian ini diedit oleh editor Jurnal Matematika UNAND \begin{abstract} \begin{center} Diterima ..... \quad Direvisi ..... \quad Dipublikasikan ..... %tanggal-tanggal tersebut \textbf{dikosongkan} saja \end{center} \bigskip %\textbf{Abstrak}. %Dalam bahasa Indonesia %Abstrak tidak melebihi 250 kata. Tuliskan intisari dari makalah anda. \bigskip \textbf{Abstract}. % Dalam bahasa Inggris \textit{The classical solution and the strong solution of a partial differential equation problem are continuously differentiable solutions. This solution has a derivative for a continuous infinity level. However, not all problems of partial differential equations can be easily obtained by strong solutions. Even the existence of a solution requires in-depth investigation. The variational formulation method can qualitatively analyze a single solution to a partial differential equation problem. This study provides an alternative method in analyzing the problem model of partial differential equations analytically. In this research, we will examine the partial differential equation modelling built from fluid dynamics modelling. } \end{abstract} \keywords{Weak solution, Hilbert space, Variational method} \section{Introduction} Partial differential equations are widely used in modelling various physical, geometric and probabilistic phenomena. The partial differential equation modelling can be in the form of linear and nonlinear equations depending on the complexity of the model to be investigated. There is no known general theory of the solvency of all partial differential equations. Such a theory is implausible, given the large variety of possible models. This issue has led to the development of research on methods of solving partial differential equations problems. Some cases of partial differential equations that are interesting and continue to be developed in various studies are the Laplace equation, Helmholtz equation, Liouville equation, wave equation, Korteweg-de Vries equation \cite{paganrev}\cite{pagani}. As for the system of partial differential equations, especially in the case of fluid mechanics, the Navier-Stokes system of nonlinear equations and their linearized form, namely the Neumann-Kelvin system of equations, is still being developed to model the dynamics of fluid motion. Given the wide variety of cases that can be modelled with partial differential equations of varying degrees of difficulty, the question of what it means to "solve" a partial differential equation may be subtle. Classical and analytical solutions are perfect solutions to the problem of partial differential equations. From these solutions, it can be seen how the dynamics of the model under study are. However, not all problems of differential equations can be obtained with analytical solutions. The informal idea in building an excellent partial differential equation problem model is a model that captures many of the desired features and meets the following three criteria \cite{brezis} (a) the problem has a solution; (b) the solution is single; (c) the solution depends continuously on the data given in the problem. The last condition is very important for problems arising from physical applications. Meanwhile, observing the existence and singularity of solutions is the main task for mathematicians and classifying and characterizing solutions. Next is to define what is meant by 'solution'. The desired solution satisfies (a) – (c). However, does the solution have to be genuinely analytic or at least infinitely differentiable? This solution is the perfect condition, but it will not be easy to achieve in some cases. Perhaps it would be wiser if the solutions of partial differential equations of order k could be differentiated as much as k times continuously. Then at least all the derivatives that appear in the problem of partial differential equations will be continuous, although it is possible that certain higher derivatives will not exist. The idea of solving a partial differential equation problem has always been an exciting topic to study. In this study, variational methods will be studied to analyze the existence and singularity of a solution qualitatively. The principle of the variational approach to solving partial differential equations is to replace the equation with an equivalent formula called variational, which is obtained by integrating the equation multiplied by any function, called the test function. The solution generated by this method is a weak solution of partial differential equations. Through the Lax-Milgram Theorem, it can be shown that the solution exists and is singular. Then through the regulatory process, it will be proven that the weak solution is a strong solution to the problem of partial differential equations. This variational method is very useful in determining the existence and singularity of solutions to the boundary condition problems of partial differential equations. One of the applications of this method can be found in the article \cite{paganrev}. In that article, they analyzed weak solutions by functional analysis to prove the existence of a single solution of the linear boundary value problem that describes the motion of the equilibrium state of a half-submerged cylinder in an ideal, incompressible heavy fluid. This study provides an alternative method in analyzing the problem model of partial differential equations analytically. The partial differential equation modelling that will be analyzed in this research is built from fluid dynamics modelling. \section{Finite element method for one dimensional case} We seek to calculate an approximation of the solution $u:[0, 1] \to \mathbb{R}$ of the following problem: \begin{align} - \epsilon u''(x) + \lambda u'(x) = f(x),\\ u(0)=u(1)=0, \end{align} with $\epsilon >0$, $\lambda$ and $f$ are given so that there is a unique solution to this problem and we want to approach this solution by a continuous function, polynomial by parts.\\ \begin{enumerate} \item \textbf{Write the variational formulation of the problem}\\ We construct the variational formulation by multiplying both side of (1) by a test function $v$ satisfying the boundary conditions and then integrate over $[0, 1]$. We thus obtain \begin{align*} -\int_0^1 - \epsilon u''(x)\,v(x) \,\mathrm{d}x+ \int_0^1 \lambda u'(x)\,v(x)\mathrm{d}x=\int_0^1 f(x)\,v(x)\,\mathrm{d}x,\quad \forall v \in H^1(0, 1). \end{align*} Integrating by parts, and since v(0)=v(1)=0 then we have\\ \begin{align}\label{VF} \epsilon \int_0^1 u'(x)\,v'(x) \,\mathrm{d}x + \lambda \int_0^1 u'(x)\,v(x)\mathrm{d}x=\int_0^1 f(x)\,v(x)\,\mathrm{d}x,\quad \forall v \in H^1(0, 1). \end{align} We obtain the bilinear and linear form: \begin{align*} B[u,v]&=\epsilon \int_0^1 u'(x)\,v'(x) \,\mathrm{d}x + \lambda \int_0^1 u'(x)\,v(x)\mathrm{d}x\\ L(v) &= \int_0^1 f(x)\,v(x)\,\mathrm{d}x \end{align*} Therefore we can express the variational formulation as follows: find $u \in H^1_0([0, 1])$ such that for all $v \in H^1_0([0, 1])$, we have $B[u,v]=L(v)$. \item \textbf{Calculation of the exact solution}\\ We suppose the function $f$ is constant not null.\\ \begin{enumerate} \item Find the exact solution of the problem.\\ $$- \epsilon u''(x) + \lambda u'(x) = f(x) \Leftrightarrow u''(x) + \frac{1}{- \epsilon}\lambda u'(x) = \frac{1}{- \epsilon} f(x)$$ - Subtitute : $$r'(x)+\frac{\lambda}{- \epsilon}r = \frac{f}{- \epsilon}$$\\ - Integrating factor:$$\mu(x) = e^{\int \frac{\lambda}{- \epsilon}dx} = e^{- \frac{\lambda}{\epsilon}x + C_1}=C_2e^{-\frac{\lambda}{\epsilon}x},\quad C_2=e^{C_1}.$$\\ - Solve for $r$: \begin{align*} r(x)&=\frac{\int C_2 e^{\frac{-\lambda}{\epsilon}x}\, \frac{f}{-\epsilon} \, dx}{{C_2 e^{\frac{-\lambda}{\epsilon}x}}}\\ &= \frac{f}{\lambda}+ C_4 e^{\frac{\lambda}{\epsilon}x},\quad C_4=C_3/C_2. \end{align*} - Solution : $$u(x)=\int r(x) \, dx=\int \frac{f}{\lambda}+ C_4 e^{\frac{\lambda}{\epsilon}x}\, dx = \frac{f}{\lambda}x+ C_4 \frac{\epsilon}{\lambda} e^{\frac{\lambda}{\epsilon}x} + C_5.$$ - Subtitute the boundary conditions, we obtain the exact solution of (1) and (2): \begin{align} u(x) = \frac{f}{\lambda}x+ \frac{f}{\lambda (1-e^{\frac{\lambda}{\epsilon}})} \left(e^{\frac{\lambda}{\epsilon}x}-1 \right). \end{align} \item Write a function to calculate the solution $u$ on a grid of points given. Draw $u$ for $f=1, \lambda \in {-1,1}$, and $\epsilon \in \{1,0.5,0.1,0.01\}.$\\ Let us see the graph of $u$ : \begin{figure}[h] \centering \includegraphics[scale=0.5]{f1l1ep1.png} \caption{$f=1, \lambda = 1, \epsilon = 1:$} \includegraphics[scale=0.5]{f1l1ep01.png} \caption{$f=1, \lambda = 1, \epsilon = 0.1:$} \end{figure} \begin{figure}[h] \centering \includegraphics[scale=0.5]{f1l1ep001.png} \caption{$f=1, \lambda = 1, \epsilon = 0.01:$} \includegraphics[scale=0.5]{f1l-1ep001.png} \caption{$f=1, \lambda = -1, \epsilon = 0.1:$} \label{} \end{figure} \end{enumerate} \newpage \item \textbf{We seek an approximate solution (finite element method P1)} \begin{enumerate} \item Write the approximate variational formulation.\\ Let us generate a mesh, let a uniform Cartesian mesh $x_i=ih, i=0,1,2,\cdots,n$ where $h=1/n$, defining the intervals $(x_{i-1},x_i)$, $i=1,2,\cdots,n$. Then, we construct a set of basis functions bassed on the mesh, such as the piecewise linear functions $i=1,2,\cdots,n-1$ \begin{align} \phi_i(x)= \left\{ \begin{array}{ll} 0 & x< x_{i-1}\\ \frac{x-x_{i-1}}{h} & x_{i-1}\leq x\leq x_{i} \\ \frac{x_{i+1}-x}{h} & x_{i}\leq x\leq x_{i+1} \\ 0 & x> x_{i+1} \\ \end{array} \right.. \end{align} \begin{figure}[h] \centering \includegraphics[scale=0.7]{P1n=20.png} \caption{basis function for $n=20$} \end{figure} Represent the approximate (FE) solution by a linear combination of the basis functions \begin{align}\label{uh} u_h(x) = \sum_{j=1}^{n-1} c_j \phi_j(x) \end{align} where the coefficients $c_j$ are the unknowns to be determined. Subtituting the approximate solution $u_h(x)$ in the variational formulation: \begin{align}\label{VF1} \epsilon \int_0^1 u_h'(x)\,v'(x) \,\mathrm{d}x + \lambda \int_0^1 u_h'(x)\,v(x)\mathrm{d}x=\int_0^1 f(x)\,v(x)\,\mathrm{d}x \end{align} Let us define a finite dimensional space on the mesh. Let the solution be in the space $V$, which is $H_0^1(0, 1)$ in the model problem. Based on the mesh, we wish to construct a subspace $V_h$ (a finite dimensional space) $\subset V$ (the solution space). $$V_h = \{v_h (x)| v_h (x) \text{is continuous piecewise linear},\, v_h(0) = v_h (1) = 0\}.$$ \item Show that the approximate solution $u_h$ can be found as the solution of a linear system $$A_h u_h=b_h.$$ Deduce that $A_h = \epsilon B_h + \lambda C_h$, where $B_h$ is a symmetric tridiagonal matrix and $C_h$ is an antisymmetric tridiagonal matrix.\\ From (\ref{VF1}), we have \begin{align}\label{VF1} \epsilon \int_0^1 \sum_{j=1}^{n-1} c_j \phi_j'(x)\,v'(x) \,\mathrm{d}x + \lambda \int_0^1 \sum_{j=1}^{n-1} c_j \phi_j'(x)\,v(x)\mathrm{d}x=\int_0^1 f(x)\,v(x)\,\mathrm{d}x.\\ \Rightarrow \epsilon\sum_{j=1}^{n-1} c_j \int_0^1 \phi_j'(x)\,v'(x) \,\mathrm{d}x + \lambda \sum_{j=1}^{n-1} c_j \int_0^1 \phi_j'(x)\,v(x)\mathrm{d}x=\int_0^1 f(x)\,v(x)\,\mathrm{d}x. \end{align} We choose the test function $v(x)$ as $\phi_1, \phi_2,\cdots,\phi_{n-1}$ successively, in the matrix vector form: \end{enumerate} \end{enumerate} \begin{small}\scriptsize \begin{align} \left( \epsilon \begin{bmatrix} b(\phi_1,\phi_1) & b(\phi_1,\phi_2) & \dots & b(\phi_1,\phi_{n-1}) \\ b(\phi_2,\phi_1) & b(\phi_2,\phi_2) & \dots & b(\phi_2,\phi_{n-1}) \\ \vdots & \vdots & \ddots & \vdots \\ b(\phi_{n-1},\phi_1) & b(\phi_{n-1},\phi_2) & \dots & b(\phi_{n-1},\phi_{n-1}) \end{bmatrix} +\lambda \begin{bmatrix} c(\phi_1,\phi_1) & c(\phi_1,\phi_2) & \dots & c(\phi_1,\phi_{n-1}) \\ c(\phi_2,\phi_1) & c(\phi_2,\phi_2) & \dots & c(\phi_2,\phi_{n-1}) \\ \vdots & \vdots & \ddots & \vdots \\ c(\phi_{n-1},\phi_1) & c(\phi_{n-1},\phi_2) & \dots & c(\phi_{n-1},\phi_{n-1}) \end{bmatrix} \right) \begin{bmatrix} c_1 \\ c_2 \\ \vdots \\ c_3 \end{bmatrix} = \begin{bmatrix} (f,\phi_1) \\ (f,\phi_2) \\ \vdots \\ (f,\phi_{n-1}) \end{bmatrix} \end{align} \end{small} where $b(\phi_i,\phi_j)=\int_0^1\phi_i'(x)\,\phi_j'(x) \,\mathrm{d}x$, $ c(\phi_i,\phi_j)=\int_0^1\phi_i'(x)\,\phi_j(x) \,\mathrm{d}x$, $(f,\phi_i)=\int_0^1f\,\phi_i(x) \,\mathrm{d}x$.\\ Then we can write the approximation solution $u_h$ can be obtained by solve the following linear system:$$A_h\,u_h=b_h$$ where $A_h=\epsilon B_h+ \lambda C_h.$\\ By standard calculation of integral, we will get $B_h$ is a symmetric tridiagonal matrix and $C_h$ is an antisymmetric tridiagonal matrix. Let us see the case: $\epsilon = 0.1$, $\lambda = 1$, $f=1$, $n=10$: \begin{align*} B_h= \left( \begin{array}{ccccccccc} 20 & -10 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ -10 & 20 & -10 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & -10 & 20 & -10 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & -10 & 20 & -10 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & -10 & 20 & -10 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & -10 & 20 & -10 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & -10 & 20 & -10 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & -10 & 20 & -10 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & -10 & 20 \\ \end{array} \right) \end{align*}, \begin{align*} C_h= \left( \begin{array}{ccccccccc} 0 & 0.5 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ -0.5 & 0 & 0.5 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & -0.5 & 0 & 0.5 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & -0.5 & 0 & 0.5 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & -0.5 & 0 & 0.5 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & -0.5 & 0 & 0.5 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & -0.5 & 0 & 0.5 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & -0.5 & 0 & 0.5 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & -0.5 & 0 \\ \end{array} \right). \end{align*} Then we have \begin{align*} A_h = \left( \begin{array}{ccccccccc} 2. & -0.5 & 0. & 0. & 0. & 0. & 0. & 0. & 0. \\ -1.5 & 2. & -0.5 & 0. & 0. & 0. & 0. & 0. & 0. \\ 0. & -1.5 & 2. & -0.5 & 0. & 0. & 0. & 0. & 0. \\ 0. & 0. & -1.5 & 2. & -0.5 & 0. & 0. & 0. & 0. \\ 0. & 0. & 0. & -1.5 & 2. & -0.5 & 0. & 0. & 0. \\ 0. & 0. & 0. & 0. & -1.5 & 2. & -0.5 & 0. & 0. \\ 0. & 0. & 0. & 0. & 0. & -1.5 & 2. & -0.5 & 0. \\ 0. & 0. & 0. & 0. & 0. & 0. & -1.5 & 2. & -0.5 \\ 0. & 0. & 0. & 0. & 0. & 0. & 0. & -1.5 & 2. \\ \end{array} \right) \end{align*} We also have \begin{align*} (f,\phi_i)= \left( \begin{array}{c} \frac{1}{10} \\ \frac{1}{10} \\ \frac{1}{10} \\ \frac{1}{10} \\ \frac{1}{10} \\ \frac{1}{10} \\ \frac{1}{10} \\ \frac{1}{10} \\ \frac{1}{10} \\ \end{array} \right) \end{align*} Hence by the calculation of $A_h^{-1}.(f,\phi_i)$, we obtain \begin{align}\label{c} \left( \begin{array}{c} c_1 \\ c_2 \\ c_3 \\ c_4 \\ c_5 \\ c_6 \\ c_7 \\ c_8 \\ c_9\\ \end{array} \right) = \left( \begin{array}{c} 0.0999661 \\ 0.199865 \\ 0.29956 \\ 0.398645 \\ 0.495902 \\ 0.587671 \\ 0.662979 \\ 0.688904 \\ 0.566678 \\ \end{array} \right). \end{align} Subtituting (\ref{c}) into (\ref{uh}) then we have the approximation solution. Let us see the table below: \begin{align} \left( \begin{array}{cccc} x & u(x) & u_h(x) & u(x)-u_h(x) \\ 0 & 0. & 0. & 0. \\ 0.1 & 0.099922 & 0.0999661 & -0.0000441427 \\ 0.2 & 0.19971 & 0.199865 & -0.000154593 \\ 0.3 & 0.299133 & 0.29956 & -0.000426202 \\ 0.4 & 0.397567 & 0.398645 & -0.00107863 \\ 0.5 & 0.493307 & 0.495902 & -0.00259449 \\ 0.6 & 0.581729 & 0.587671 & -0.00594212 \\ 0.7 & 0.650256 & 0.662979 & -0.0127232 \\ 0.8 & 0.664704 & 0.688904 & -0.0242 \\ 0.9 & 0.532149 & 0.566678 & -0.0345287 \\ 1 & 0. & 0. & 0. \\ \end{array} \right) \end{align} For $n=20$, we can obtain the solution with the similar calculation as follows: \begin{align} \left( \begin{array}{ccc} u(x) & u_h(x) & u(x)-u_h(x) \\ 0. & 0. & 0. \\ 0.0499705 & 0.0499756 & -5.0779565\times10^{-6}\\ 0.099922 & 0.099935 & -0.0000130127 \\ 0.149842 & 0.149867 & -0.0000253658 \\ 0.19971 & 0.199754 & -0.0000445175 \\ 0.249492 & 0.249566 & -0.0000740682 \\ 0.299133 & 0.299253 & -0.000119414 \\ 0.348542 & 0.34873 & -0.000188551 \\ 0.397567 & 0.39786 & -0.000293164 \\ 0.445958 & 0.446408 & -0.000450015 \\ 0.493307 & 0.49399 & -0.000682575 \\ 0.538936 & 0.539958 & -0.0010226 \\ 0.581729 & 0.58324 & -0.00151086 \\ 0.619847 & 0.622042 & -0.00219529 \\ 0.650256 & 0.653379 & -0.00312278x \\ 0.667957 & 0.672274 & -0.00431704 \\ 0.664704 & 0.670432 & -0.00572785 \\ 0.626905 & 0.634029 & -0.00712355 \\ 0.532149 & 0.540023 & -0.00787414 \\ 0.343487 & 0.350015 & -0.00652742 \\ 0. & 0. & 0. \\ \end{array} \right) \end{align} \begin{figure}[h] \centering \includegraphics[scale=0.7]{graph1.png} \caption{Comparison result between exact solution and approximation solution for $P1(n=10,n=20)$} \end{figure} Study of the error. Fix $\epsilon = 0.1, \lambda = 1, f=1$. For $n$ ranging from 10 to 100 (in steps of 10), plot the curves $\log n \to \log ||e_n||_\infty$ where $e_n$ is the error vector defined by $e_n (k) = u_h(k) - u(x_k).$ Deduce a decay law of $||e_n||_\infty$ of the form $||e_n||_\infty\sim Constant/n^{s_1}$ where $s_1$ is to be determined. \textbf{4. We seek an approximate solution (finite element method P2)} The basis function for P2 is quadratic function. Let us generate a mesh, let a uniform Cartesian mesh $x_i=ih, i=0,1,2,\cdots,n$ where $h=1/n$, defining the intervals $(x_{i-1},x_i)$, $i=1,2,\cdots,n$. Then, we construct two sets of functions bassed on the mesh, such as the quadratic functions $i=1,2,\cdots,n-1$. The first function is: \begin{align} p_i(x)= \left\{ \begin{array}{ll} 0 & x< x_{i-1}\\ \frac{(x-x_{i-0.5})(x-x_{i-1})}{0.5h^2} & x_{i-1}\leq x\leq x_{i} \\ \frac{(x-x_{i+0.5})(x-x_{i+1})}{0.5h^2} & x_{i}\leq x\leq x_{i+1} \\ 0 & x> x_{i+1} \\ \end{array} \right.. \end{align} The second function is for fill the middle value of the first function \begin{align} q_j(x)= \left\{ \begin{array}{ll} -\frac{(x-x_{i})(x-x_{i+1})}{0.25h^2} & x_{i}\leq x\leq x_{i+1} \\ 0 & \text{otherwise} \\ \end{array} \right. \end{align} for $j=0,1,2,...,n-1$ Then we construct the basis function as the combination of that two functions such as: $$\psi_j=\{\psi_1=p_1, \psi_2=p_2, \cdots, \psi_{n-1}=p_{n-1}, \psi_{n}=q_0, \psi_{n+1}=q_1,\cdots,\psi_{2n-1}=q_{n-1}\}.$$ We can plot the basis function as follows (for n=10): \begin{figure}[h] \centering \includegraphics[scale=0.5]{figp2.png} \caption{basis function for $P2$ case} \end{figure} Represent the approximate (FE) solution by a linear combination of the basis functions \begin{align}\label{uh} u_h(x) = \sum_{j=1}^{2n-1} c_j \psi_j(x) \end{align} where the coefficients $c_j$ are the unknowns to be determined. For the case, $n=4,\epsilon = 0.1, \lambda=1,f=1$ we obtain these results: \begin{align*} B_h=\left( \begin{array}{ccccccc} 18.667 & 1.333 & 0 & -10.667 & -10.667 & 0 & 0 \\ 1.333 & 18.667 & 1.333 & 0 & -10.667 & -10.667 & 0 \\ 0 & 1.333 & 18.667 & 0 & 0 & -10.667 & -10.667 \\ -10.667 & 0 & 0 & 21.333 & 0 & 0 & 0 \\ -10.667 & -10.667 & 0 & 0 & 21.333 & 0 & 0 \\ 0 & -10.667 & -10.667 & 0 & 0 & 21.333 & 0 \\ 0 & 0 & -10.667 & 0 & 0 & 0 & 21.333 \\ \end{array} \right). \end{align*} \begin{align*} C_h=\left( \begin{array}{ccccccc} 0 & -0.167 & 0 & -0.667 & 0.667 & 0 & 0 \\ 0.167 & 0 & -0.167 & 0 & -0.667 & 0.667 & 0 \\ 0 & 0.167 & 0 & 0 & 0 & -0.667 & 0.667 \\ 0.667 & 0 & 0 & 0 & 0 & 0 & 0 \\ -0.667 & 0.667 & 0 & 0 & 0 & 0 & 0 \\ 0 & -0.667 & 0.667 & 0 & 0 & 0 & 0 \\ 0 & 0 & -0.667 & 0 & 0 & 0 & 0 \\ \end{array} \right). \end{align*} \begin{align*} A_h=0.1*B_h+C_h=\left( \begin{array}{ccccccc} 1.8667 & -0.0337 & 0. & -1.7337 & -0.3997 & 0. & 0. \\ 0.3003 & 1.8667 & -0.0337 & 0. & -1.7337 & -0.3997 & 0. \\ 0. & 0.3003 & 1.8667 & 0. & 0. & -1.7337 & -0.3997 \\ -0.3997 & 0. & 0. & 2.1333 & 0. & 0. & 0. \\ -1.7337 & -0.3997 & 0. & 0. & 2.1333 & 0. & 0. \\ 0. & -1.7337 & -0.3997 & 0. & 0. & 2.1333 & 0. \\ 0. & 0. & -1.7337 & 0. & 0. & 0. & 2.1333 \\ \end{array} \right) \end{align*} \begin{align*} A_h^{-1}\left( \begin{array}{ccccccc} 0.901412 & 0.0875424 & 0.00781085 & 0.732563 & 0.240035 & 0.0227499 & 0.00146346 \\ 0.893776 & 0.981158 & 0.0875424 & 0.726358 & 0.964832 & 0.254976 & 0.0164022 \\ 0.814177 & 0.893776 & 0.901412 & 0.661669 & 0.878904 & 0.900023 & 0.168891 \\ 0.168891 & 0.0164022 & 0.00146346 & 0.606012 & 0.0449735 & 0.00426248 & 0.000274197 \\ 0.900023 & 0.254976 & 0.0227499 & 0.731435 & 0.844603 & 0.0662615 & 0.00426248 \\ 0.878904 & 0.964832 & 0.240035 & 0.714272 & 0.948777 & 0.844603 & 0.0449735 \\ 0.661669 & 0.726358 & 0.732563 & 0.537729 & 0.714272 & 0.731435 & 0.606012 \\ \end{array} \right). \end{align*} \begin{align*} F= \left( \begin{array}{c} 0.083 \\ 0.083 \\ 0.083 \\ 0.167 \\ 0.167 \\ 0.167 \\ 0.167 \\ \end{array} \right) \end{align*} \begin{align*} c=A_h^{-1}.F=\left( \begin{array}{c} 0.249199 \\ 0.490635 \\ 0.652362 \\ 0.124973 \\ 0.372729 \\ 0.599242 \\ 0.608447 \\ \end{array} \right) \end{align*} Then we have \begin{align}\label{uh} u_h(x) = \sum_{j=1}^{7} c_j \psi_j(x). \end{align} \begin{align*} \left( \begin{array}{cccc} x &u(x) & u_h(x) & u(x)-u_h(x)\\ 0 &0. & 0 & 0 \\ 1/8 &0.124887 & 0.124973 & -0.0000860471 \\ 2/8& 0.249492 & 0.249199 & 0.000293245 \\ 3/8&0.373115 & 0.372729 & 0.000385755 \\ 4/8&0.493307 & 0.490635 & 0.00267265 \\ 5/8&0.601527 & 0.599242 & 0.00228501 \\ 6/8&0.667957 & 0.652362 & 0.0155949 \\ 7/8&0.588528 & 0.608447 & -0.0199193 \\ 1&0 & 0 & 0 \\ \end{array} \right) \end{align*} \begin{figure} \centering \includegraphics[scale=0.5]{P2sim.png} \caption{Comparison result between exact solution and approximation solution for $P2(n=4)$} \end{figure} Now let us check for $n=10$: \begin{align*} \left( \begin{array}{cccc} x&u(x) & u_h(x) & u(x)-u_h(x) \\ 0&0. & 0 & 0 \\ 0.05&0.0499705 & 0.0499826 & -0.0000120677 \\ 0.1&0.099922 & 0.0999511 & -0.000029145 \\ 0.15&0.149842 & 0.149892 & -0.0000502498 \\ 0.2&0.19971 & 0.199786 & -0.0000764664 \\ 0.25&0.249492 & 0.249599 & -0.000106273 \\ 0.3&0.299133 & 0.299273 & -0.000139483 \\ 0.35&0.348542 & 0.348719 & -0.000177509 \\ 0.4&0.397567 & 0.397779 & -0.000212591 \\ 0.45&0.445958 & 0.446217 & -0.000258685 \\ 0.5&0.493307 & 0.493591 & -0.000283661 \\ 0.55&0.538936 & 0.539276 & -0.000340579 \\ 0.6&0.581729 & 0.582057 & -0.00032816 \\ 0.65&0.619847 & 0.620259 & -0.000411977 \\ 0.7&0.650256 & 0.650558 & -0.000302309 \\ 0.75&0.667957 & 0.668437 & -0.000480028 \\ 0.8&0.664704 & 0.664853 & -0.000149252 \\ 0.85&0.626905 & 0.627571 & -0.000666252 \\ 0.9&0.532149 & 0.532032 & 0.000116878 \\ 0.95&0.343487 & 0.345065 & -0.00157802 \\ 1&0. & 0 & 0 \\ \end{array} \right) \end{align*} \begin{figure}[h] \centering \includegraphics[scale=0.6]{graph2.png} \caption{Comparison result between exact solution and approximation solution for $P2(n=10)$} \end{figure} \section{Variational method for two dimensional case} Let $\Omega$ an open bounded related of $\mathbb{R}^N$, $N \geq 1$. For $v \in L^2(\Omega)$, we denote \begin{align} \bar{v}=\frac{1}{|\Omega|}\int_{\Omega} v(x) \mathrm{d}x \end{align} is the average of $v$ in $\Omega$ and $|\Omega|$ is the measure of $\Omega$. We consider the Hilbert space $V$, \begin{align} V= \{v \in H^1(\Omega), \bar{v} =0 \} \end{align} endowed with norm: \begin{align} \lVert v \rVert^2=\int_{\Omega} |\nabla v|^2 \mathrm{d}x \end{align} \begin{theorem} There exists a unique function $u \in V$ such that \begin{align} \label{VP1} \int_{\Omega} \nabla u \nabla v \mathrm{d}x = \int_{\Omega} gv \mathrm{d}x \quad \forall v\in V \end{align} for $g \in L^2(\Omega)$ \end{theorem} \begin{proof} We will apply the Lax-Milgram theorem \cite{evans} by checking the continuity of bilinear and linear form and also coercivity of the bilinear form. From \eqref{VP1} we have the bilinear form is \begin{align} B[u,v] = \int_{\Omega} \nabla u \nabla v \mathrm{d}x \quad \forall u, v\in V \end{align} and the linear form is \begin{align} L(v) = \int_{\Omega} gv \mathrm{d}x \quad \forall u, v\in V\, \text{and} \, g \in L^2(\Omega) \end{align} First we will check the continuity of $B$, \begin{align} \left|B[u,v]\right| &= \left| \int_{\Omega} \nabla u \nabla v \mathrm{d}x \right|\\ &= \lVert \nabla u \nabla v \rVert_{L^1(\Omega)}\\ & \leq \lVert \nabla u \rVert_{L^2(\Omega)}\, \lVert \nabla v \rVert_{L^2(\Omega)}\\ &= \lVert u \rVert_V\, \lVert v \rVert_V \end{align} We obtain $\left|B[u,v]\right| \leq \lVert u \rVert_V\, \lVert v \rVert_V$, $B$ continue in $V$.\\ \\* Next we will check the coercivity of $B$ \begin{align} B[u,u] = \int_{\Omega} |\nabla u|^2 \mathrm{d}x = \lVert v \rVert_V^2 \end{align} The coercivity is satistied by the norm of $V$.\\ \\* Finally we obtain the continuity of $L$ by Cauchy-Schwarz inequality, \begin{align} |L(v)| &= \left|\int_{\Omega} gv \mathrm{d}x \right|\\ &\leq \lVert g \rVert_{L^2(\Omega)}\, \lVert v \rVert_{L^2(\Omega)}\\ &\leq \lVert g \rVert_{L^2(\Omega)}\, \lVert v \rVert_V \end{align} According to the Lax-Milgram theorem, we obtain that there exists a unique function $u\in V$ such that for all $v \in V$, we have $$B[u,v] = L(v)$$ \end{proof} %\textbf{Q2}. Show that the unique solution of the problem \eqref{VP1} %verify: %\begin{align} % \int_{\Omega} \nabla u \nabla v \mathrm{d}x = \int_{\Omega} (g-\bar g)v \mathrm{d}x \quad \forall v\in H^1(\Omega) % \end{align} %Deduce for $u \in H^2(\Omega)$ and if $g \in C^{1,\alpha}(\bar \Omega)$ then $ u \in C^3(\bar \Omega)$. %\\* %\textbf{Solution}\\ %Since $v \in H^1(\Omega)$ and $\bar{v} = \frac{1}{|\Omega|}\int_{\Omega} v(x) \mathrm{d}x$ is the average %of $v$ then $v-\bar{v} \in V$, thereby we deduce, %\begin{align} % \int_{\Omega} \nabla u \nabla (v-\bar{v}) \mathrm{d}x &= \int_{\Omega} g(v-\bar{v}) \mathrm{d}x\\ %\int_{\Omega} \nabla u \nabla v\mathrm{d}x-\int_{\Omega} \nabla u \nabla \bar{v} \mathrm{d}x &= \int_{\Omega} gv\mathrm{d}x-\int_{\Omega} g\bar{v} \mathrm{d}x % \end{align} %Let us have a look that $\bar{v} = \frac{1}{|\Omega|}\int_{\Omega} v(x) \mathrm{d}x \in \mathbb{R}$ is a constant,$\nabla \bar{v}=0$. Then we have %\begin{align} %\int_{\Omega} \nabla u \nabla v\mathrm{d}x &= \int_{\Omega} gv\mathrm{d}x-\bar{v}\int_{\Omega} g \mathrm{d}x\\ %&= \int_{\Omega} gv\mathrm{d}x-\bar{v}|\Omega| \bar{g} \\ %&= \int_{\Omega} gv\mathrm{d}x-\bar{g}\int_{\Omega} v\mathrm{d}x \\ %&= \int_{\Omega} (g-\bar{g}) v\mathrm{d}x % \end{align} % %Next let us see Theorem 5.2 in Rakotoson,\\ %\textbf{Theorem 5.2} %\textit{Let $\Omega$ is an open bounded set of class $C^{m+2}, m \geq 0, f \in C^{m,a}(\bar{\Omega})$ of the problem: %\begin{align} % \int_{\Omega} \nabla u \nabla v \mathrm{d}x +\int_{\Omega} u v \mathrm{d}x= \int_{\Omega} fv \mathrm{d}x \quad \forall v\in H^1(\Omega) % \end{align} %then $u \in C^{m+2,a}(\bar{\Omega})$ % } %\\* %From (21), if we define $F=g-\bar{g}+u \in L^2(\Omega)$ then we can rewrite the equation as follows: %\begin{align} % \int_{\Omega} \nabla u \nabla v \mathrm{d}x +\int_{\Omega} u v \mathrm{d}x= \int_{\Omega} Fv \mathrm{d}x \quad \forall v\in H^1(\Omega) % \end{align} % Let us see the equations we have, % \begin{align*} % \int_{\Omega} \nabla u \nabla v \mathrm{d}x &= \int_{\Omega} gv \mathrm{d}x \quad \forall v\in V\\ % \int_{\Omega} \nabla u \nabla v \mathrm{d}x &= \int_{\Omega} (g-\bar{g})v \mathrm{d}x \quad \forall v\in H^1(\Omega) % \end{align*} %From those, for all $v \in C_0^{\infty}(\Omega) \subset H^1(\Omega)$ we have % \begin{align*} % \int_{\Omega} \nabla u \nabla v \mathrm{d}x &= \int_{\Omega} (g-\bar{g})v \mathrm{d}x % \end{align*} %From those we have, % \begin{align} % -\Delta u = g-\bar{g} \, \text{(in the distribution sense)}\ % \end{align} %with $g-\bar{g} \in L^2(\Omega)$ and $u \in V$, we obtain $u \in H^2(\Omega)$. Also, by the Sobolev embedding theorems, %$H^2(\Omega)\subset C^0(\bar{\Omega})$ then $u \in C^0(\bar{\Omega}), F \in C^0(\bar{\Omega})$ by Theorem 5.2, it implies $u \in C^2(\bar{\Omega})$. %Futhermore, $F=(g-\bar{g})+u$ with $(g-\bar{g}) \in C^1(\bar{\Omega})$ and $u \in C^2(\bar{\Omega})$ then we have, %$F \in C^1(\bar{\Omega})$, and we apply the Theorem 5.2 for second time then we have $ u \in C^3(\bar{\Omega})$ %\\ %\\* %\textbf{Q3}. Let $g\in C^{1,\alpha}(\bar{\Omega}) (0<\alpha<1)$ with zero average (i.e $\bar{g}=0)$. %Show that there exists at least one vector field $\vec{U}$ of class $C^2$ in $\Omega$ such that: %\begin{align} % \mathrm{div}(\vec{u}) &= g \quad \text{in}\; \Omega\\ % \vec{u} \cdot \vec{n} &= 0 \quad \text{on} \;\partial\Omega %\end{align} %where $\vec{n}(x)$ is the exterior normal of point $x$ from the edge of $\partial\Omega$. %\\* %\textbf{Solution}\\ %From Q1 and Q2 we have if $g \in C^{1,\alpha}(\bar{\Omega})$ then there exists $u \in C^3(\bar{\Omega})$ verify % \begin{align}\label{24.3} % \int_{\Omega} \nabla u \nabla v \mathrm{d}x = \int_{\Omega} gv \mathrm{d}x \quad \forall v\in H^1(\Omega) % \end{align} %If $v \in C_0^{\infty}(\Omega)$ then we have %\begin{align} % \big \langle \nabla u, \nabla v \big \rangle = \big \langle g, v \big \rangle\\ % \big \langle -\Delta u, v \big \rangle = \big \langle g, v \big \rangle % \end{align} %We obtain $-\Delta u = g$ in the distribution sense. Thus if we take a test function $v\in H^1(\Omega)$ %and multiply by $-\Delta u = g$ and then take the integration over $\Omega$, we have %\begin{align} % -\int_{\Omega}\Delta u v \mathrm{d}\Omega = \int_{\Omega}g v \mathrm{d}\Omega % \end{align} %Integrating by parts, we obtain %\begin{align}\label{24.4} % -\int_{\partial\Omega}\nabla u \cdot n \, v \mathrm{d}(\partial\Omega)+\int_{\Omega}\nabla u \,\nabla v \mathrm{d}(\Omega) = \int_{\Omega}g v \mathrm{d}\Omega % \end{align} %where $n$ is the exterior normal of point $x$ from the edge of $\partial\Omega$.\\ %\\* %If we comparing \eqref{24.4} with \eqref{24.3}, we have: %\begin{align} %-\int_{\partial\Omega}\nabla u \cdot n \, v \mathrm{d}(\partial\Omega)=0 \quad \forall v \in H^1(\Omega) %\end{align} %Let us see this theorem,\\ %\textbf{Theorem} \textit{Let $f \in L^2(\partial\Omega)$, if %\begin{align} % \int_{\partial \Omega} fg=0, \; \forall g\in L^2(\partial\Omega)\; then\; f=0 \; in \;L^2(\partial\Omega) %\end{align} %} %\\* According to that theorem, since $trace \, v\in H^{1/2}(\partial\Omega)$ and $H^{1/2}(\partial\Omega)$ is dense in $L^2(\partial \Omega)$ then we obtain %$\nabla u \cdot n =0$ in $L^2(\partial\Omega)$.\\ % %Futhermore, for $u \in C^3(\bar{\Omega})$ we define the vector field $\vec{U}=-\nabla u$ in $\Omega$ such that these following conditions are satisfied: %\begin{align} % \mathrm{div}(\vec{U}) &= g \quad \text{in}\; \Omega \label{vec1}\\ % \vec{U} \cdot \vec{n} &= 0 \quad \text{on} \;\partial\Omega \label{vec2} %\end{align} %\\ %\\* %\textbf{Q4}. Consider this following problem: Find $(\vec{u},p)\in H^1(\Omega) \times L^2(\Omega)/\mathbb{R}$ such that: %\begin{align} % -\Delta \vec{u} + \nabla p &= \vec{f} \quad \in L^2(\Omega)^N\\ % \vec{u} \cdot \vec{n} &= 0\quad \text{on}\, \partial \Omega\\ % \mathrm{div}(\vec{u}) &= g \quad \text{where} \, g \in C^{1,\alpha}(\bar{\Omega}), \bar{g}=0 %\end{align} %by using question 3, find the variational formulation of the system above and solve it. %\\ %\\* %\textbf{Solution}\\ %Let $\vec{U}$ is a vector of class $C^2$ and satisfied \eqref{vec1} and \eqref{vec2}. Then we pose %$\vec{w}=\vec{u}-\vec{U}$, we will check what problem is satisfied by $\vec{w}$. %\begin{align} % \vec{w} \cdot \vec{n} = (\vec{u}-\vec{U})\cdot \vec{n}=\vec{u}\cdot \vec{n} -\vec{U}\cdot \vec{n}= 0\\ % \mathrm{div}(\vec{w}) = \mathrm{div}(\vec{u}-\vec{U})=\mathrm{div}(\vec{u})-\mathrm{div}(\vec{U})=g-g=0 %\end{align} %Let us see equation (36), %\begin{align} % -\Delta \vec{u} + \nabla p &= \vec{f} \\ % -\Delta (\vec{w}+\vec{U}) + \nabla p &= \vec{f} \\ %-\Delta \vec{w}-\Delta \vec{U} + \nabla p &= \vec{f} % \end{align} %since $\vec{U} \in C^2(\bar{\Omega})$ then $\Delta\vec{U} \in C(\bar{\Omega})$, we have %\begin{align} % -\Delta \vec{w} + \nabla p &= \vec{f}+\Delta \vec{U} %\end{align} %For $\vec{F_1} = \vec{f}+\Delta \vec{U} \in L^2(\Omega)^N$, we obtain the problem which satisfied %by $\vec{w}$. Find $(\vec{w},p)\in H^1(\Omega) \times L^2(\Omega)/\mathbb{R}$ such that: %\begin{align} % -\Delta \vec{w} + \nabla p &= \vec{F_1} \quad \in L^2(\Omega)^N\\ % \vec{w} \cdot \vec{n} &= 0\quad \text{on}\, \partial \Omega\\ % \mathrm{div}(\vec{w}) &= 0 \quad \text{in}\, \Omega %\end{align} %We define %\begin{align} % V=\left\{\vec{v}\in \left(H^1(\Omega)\right)^N \;|\;\mathrm{div}(\vec{v}) = 0 \; \text{in}\,\Omega % \;\text{and}\;\vec{v} \cdot \vec{n} = 0\; \text{on}\; \partial \Omega\right\} % \end{align} %and endowed with norm %\begin{align} % \lVert v \rVert_V^2=\int_{\Omega} |v|^2 \, \mathrm{d}\Omega + \int_{\Omega} |\nabla v|^2 \, \mathrm{d}\Omega = \lVert v \rVert_{H^1(\Omega)}^2 %\end{align} %We will show that $V$ is a Hilbert space. Let a mapping %\begin{align} % F:(H^1)^N\,&\longrightarrow \, L^2(\Omega) \times L^2(\partial\Omega)\\ % v \, &\longrightarrow \, (\mathrm{div}(v), v \cdot n) %\end{align} %Since $\mathrm{div}(\vec{v}) = 0 \; \text{in}\,\Omega \;\text{and}\;\vec{v} \cdot \vec{n} = 0\; \text{on}\; \partial \Omega$ %then we have %\begin{align} % ker F = V = F^{-1}(\{0\}). %\end{align} %Next we will check the continuity of $F$. %\begin{align} % \lVert F(\vec{v}) \rVert_{L^2(\Omega) \times L^2(\partial\Omega)}^2 &= \lVert (\mathrm{div}(\vec{v}), \vec{v} \cdot \vec{n}) \rVert_{L^2(\Omega) \times L^2(\partial\Omega)}^2\\ %&= \lVert \mathrm{div}(\vec{v}) \rVert_{L^2(\Omega)}^2 + \lVert \vec{v} \cdot \vec{n} \rVert_{L^2(\partial\Omega)}^2 % \end{align} %We will check each term separately. %\begin{align} % \lVert \mathrm{div}(\vec{v}) \rVert_{L^2(\Omega)}^2 &= \int_{\Omega} \left(\sum_{i=1}^N \partial_i v_i \right)^2 \mathrm{d}\Omega\\ % &\leq \int_{\Omega} C_1 \sum_{i=1}^N (\partial_i v_i )^2\mathrm{d}\Omega, \quad (C \geq 3)\\ % &\leq \sum_{i=1}^N \int_{\Omega} v_i ^2 \mathrm{d}\Omega+\sum_{i=1}^N \int_{\Omega}(\partial_i v_i )^2 \mathrm{d}\Omega %\\ &= \lVert (\vec{v}) \rVert_{H^1(\Omega)}^2 %\end{align} %also %\begin{align} % \lVert \vec{v} \cdot \vec{n} \rVert_{L^2(\partial\Omega)}^2 &= \lVert \sum_{i=1}^N v_i n_i \rVert_{L^2(\partial\Omega)}^2\\ % &= \int_{\Omega} \left(\sum_{i=1}^N v_i n_i \right)^2 \mathrm{d}\Omega\\ % &\leq \int_{\Omega} C_2 \sum_{i=1}^N ( v_i n_i )^2\mathrm{d}\Omega\\ % &\leq \int_{\Omega} C_2 \sum_{i=1}^N v_i^2 \mathrm{d}\Omega \quad \text{since} \, \lVert \vec{n} \rVert =1 \, \text{then} \, n_i \leq 1, \forall i\\ % &=C_2\lVert \vec{v} \rVert_{L^2(\partial \Omega)}^2 % \\ &\leq C_3\lVert (\vec{v}) \rVert_{(H^1(\Omega))^N}^2 %\end{align} %We obtain %\begin{align} % \lVert F(\vec{v}) \rVert_{L^2(\Omega) \times L^2(\partial\Omega)}^2 \leq C_4 \lVert (\vec{v}) \rVert_{(H^1(\Omega))^2}^2 %\end{align} %which mean that $F$ continue then $ker \, F$ is closed. Since $ker \, F=V$ closed and $V \subset H^1(\Omega)$ %imply $(V, \lVert \cdot \rVert_{H^1(\Omega)})$ is a Hilbert space. % % %Moreover let us take a test function $\vec {v} \in V$ and multiply by (45) then integrate over $\Omega$ %\begin{align} % -\int_{\Omega}\Delta \vec{w}\vec {v} \mathrm{d}\Omega + \int_{\Omega}\nabla p\,\vec {v}\mathrm{d}\Omega &= \int_{\Omega}\vec{F_1}\,\vec {v}\mathrm{d}\Omega %\\-\int_{\Omega}\sum_{i=1}^N (\partial_i^2 w_i)v_i \mathrm{d}\Omega + \int_{\Omega}\sum_{i=1}^N (\partial_i p)v_i\mathrm{d}\Omega &= \int_{\Omega}\vec{F_1}\,\vec {v}\mathrm{d}\Omega % \end{align} %Integrating by parts, %\begin{align} % -\int_{\partial\Omega}\sum_{i=1}^N (\partial_i w_i)n_i \, v_i \mathrm{d}\Omega + \int_{\Omega}\sum_{i=1}^N (\partial_i w_i)\partial_i v_i \mathrm{d}\Omega % + \int_{\partial\Omega}\sum_{i=1}^N p\,v_i\,n_i\mathrm{d}\Omega - \int_{\Omega}\sum_{i=1}^N p\,\partial_i v_i\,\mathrm{d}\Omega % = \int_{\Omega}\vec{F_1}\,\vec {v}\mathrm{d}\Omega %\end{align} %since $\mathrm{div}(\vec{v}) = 0 \; \text{in}\,\Omega % \;\text{and}\;\vec{v} \cdot \vec{n} = 0\; \text{on}\; \partial \Omega$, we obtain %\begin{align} % \int_{\Omega}\sum_{i=1}^N (\partial_i w_i)(\partial_i v_i) \mathrm{d}\Omega= % \int_{\Omega}\vec{F_1}\,\vec {v}\mathrm{d}\Omega %\end{align} %Hence we have the variational formulation as follows: Find $\vec{w}\in V$ such that for all $\vec{v}\in V$, we have %\begin{align} % \int_{\Omega}\nabla \vec{w} \,\nabla\vec{v} \mathrm{d}\Omega&= % \int_{\Omega}\vec{F_1}\,\vec {v}\mathrm{d}\Omega\\ % B[\vec{w},\vec{v}]&=L(\vec{v}) %\end{align} %We will check the continuity of $B$, %\begin{align} % \left|B[\vec{w},\vec{v}]\right|&=\left|\int_{\Omega}\nabla \vec{w} \,\nabla\vec{v} \mathrm{d}\Omega\right|\\ %&\leq \lVert \nabla \vec{w} \rVert_{L^2(\Omega)}^2\,\lVert \nabla \vec{v} \rVert_{L^2(\Omega)}^2\\ %&\leq \lVert \vec{w} \rVert_V^2\,\lVert \vec{v} \rVert_V^2 %\end{align} %Next we will check the coercivity of $B$, %\begin{align} % B[\vec{w},\vec{w}] &= \int_{\Omega}|\nabla \vec{w}|^2 \mathrm{d}\Omega\\ % &=\lVert \nabla \vec{w} \rVert_{L^2(\Omega)}^2\\ % &\geq C\lVert \vec{w} \rVert_{L^2(\Omega)}^2\quad \text{Poincare inequality}\\ % &\geq \frac{C}{2} \lVert \vec{w} \rVert_{L^2(\Omega)}^2 + \frac{1}{2}\lVert \nabla \vec{w} \rVert_{L^2(\Omega)}^2\\ % &\geq \text{min}\left\{\frac{C}{2},\frac{1}{2}\right \} \lVert \vec{w} \rVert_V %\end{align} %Finally we will show the continuity of $L$, %\begin{align} % L(\vec{v})&=\int_{\Omega}\vec{F_1}\,\vec {v}\mathrm{d}\Omega\\ % &\leq \lVert \vec{F_1} \rVert_{L^2(\Omega)}^2\, \lVert \vec{v} \rVert_{L^2(\Omega)}^2\\ % &\leq \lVert \vec{F_1} \rVert_{L^2(\Omega)}^2\, \lVert \vec{v} \rVert_V^2 %\end{align} By using the FreeFem software, we can simulate the solution as follows: \begin{figure}[h] \centering \includegraphics[scale=0.3]{poisson1.png} \caption{Mesh of finite element method} \includegraphics[scale=0.3]{poisson3.png} \caption{Solution of variational formulation} \end{figure} \section{Acknowledgements} Financial support from the DIPA Fakultas MIPA Universitas Tanjungpura 2021 is gratefully acknowledged.\\ \bibliographystyle{abbrv} \bibliography{sample} %\begin{thebibliography}{0} %\bibitem{A} Nama Belakang Penulis Pertama, Singkatan Nama Depan., Nama Belakang Penulis Kedua, Singkatan Nama Depan., Tahun Penerbitan Artikel, Judul Artikel, \emph{Nama Jurnal}, \textbf{Volume} : 11 -- 22 %\bibitem{B} Nama Belakang Penulis Pertama, Singkatan Nama Depan., Nama Belakang Penulis Kedua, Singkatan Nama Depan., Tahun Penerbitan Buku, \emph{Judul Buku}, Edisi ke-, Nama Penerbit, Kota Penerbit %\bibitem{C} Nama Belakang Penulis Pertama, Singkatan Nama Depan., Tahun, \emph{Judul skripsi/tesis/disertasi}, Skripsi/Tesis/Disertasi di Nama Perguruan Tinggi, tidak diterbitkan %\end{thebibliography} \end{document}