Finite Element Method Poisson Equation in 1D with Inhomogenous Boundary Conditions 
Im trying to solve the Poisson equation in 1D: $$-u_{xx} = f(x), \hspace{6mm} u(a) = d1, \hspace{2mm} u(b) = d2$$Assuming a uniform partition such that $x_n = a + nh$, where $h = (b-a)/N$ and $n \in [0,N]$, and then discretising the problem with linear finite elements to obtain a linear equation system $\mathbf{A u} = \mathbf{f}$. I Want to find the analytical expressions for $\mathbf{A}$ and $\mathbf{f}$.

I found the general expression for $\mathbf{A}$ before incorporating boundary conditions to be $$ \mathbf{A} = \frac{1}{h}\begin{bmatrix} 1 & -1 & 0 & \ldots & \ldots & 0 \\ -1 & 2 & -1 & \ddots   & & 0 \\ 0 & -1 & 2 & \ddots & \ddots & 0 \\ \vdots & \ddots & \ddots & \ddots & \ddots  & 0 \\ \vdots &  & \ddots & -1 & 2 & -1 \\  0 & \ldots & \ldots & 0 & -1 & 1 \end{bmatrix}$$ and correspondingly for $\mathbf{f}$: $$\mathbf{f} = \big[\frac{1}{2}f(x_0), f(x_1), f(x_2), ... , \frac{1}{2}f(x_{N+1}) \big]^T$$ My trouble is with incorporating the inhomogeneous boundary conditions, I can't find any clear examples of how to do this online despite looking up a ton of sources. Is anyone able to help?
 A: Dirichlet boundary conditions are fairly easily implemented: You know your values at both ends of the interval. So you can set up "dummy" entries in your system matrix that immediately set the coefficient of the basis functions at the ends of the interval to the prescribed value:
$$ \mathbf{A} = \frac{1}{h}\begin{bmatrix} \color{red} h & \color{red} 0 & 0 & \ldots & \ldots & 0 \\ -1 & 2 & -1 & \ddots   & & 0 \\ \vdots & \ddots & \ddots & \ddots & \ddots  & 0 \\ \vdots &  & \ddots & -1 & 2 & -1 \\  0 & \ldots & \ldots & 0 & \color{red} 0 & \color{red} h \end{bmatrix}$$
and set the corresponding entries in the RHS to the known values:
$$\mathbf{f}= \begin{bmatrix} \color{red} {d1} \\ f(x_1) \\ \vdots \\ f(x_N) \\ \color{red}{d2} \end{bmatrix} $$
Then you can solve for the coefficients $u_1, \dots, u_N$ via
$$\mathbf{A} \mathbf{u} = \mathbf{f}.$$
But keep in mind that in the FEM the RHS $\mathbf{f}$ is not simply given by the point evaluation $f(x_i)$, but instead by
$$f_i = \int_{x_{i-1}}^{x_{i +1}} f(x) \phi_i(x) \mathrm d x \:, $$
i.e., by the integration of $f$ over the support of the $i$'th test function $\phi_i(x)$.
