Finding an $H^1$ function given its values on $\partial\Omega$ Background
I've met this problem when I was trying to convert a elliptic PDE problem
into  the corresponding variational problem in order to apply finite element method.
The PDE is an elliptic PDE with non-zero Dirichlet boundary condition:
Denote
$$
Lu=-\nabla\cdot(a\nabla u)+bu
$$
Then the equation is 
$$
\left\{\!\!
\begin{aligned}
&Lu=f,x\in\Omega\\
&u|_{\partial \Omega}=g
\end{aligned}
\right.
$$
When $g\equiv0$, I know the corresponding variational problem is
find $u\in H_0^1(\Omega)$, such that
$$
a(u,v)=(f,v), \forall v\in H_0^1(\Omega)
$$
where
$$
\begin{aligned}
a(u,v)&:=\int_\Omega a\nabla u\cdot\nabla v\,dx+\int_\Omega buv\,dx\qquad \\
(f,v)&:=\int_\Omega fv\,dx,\qquad \forall u,v\in H_0^1(\Omega)
\end{aligned}
$$
(This is actually the weak form of the original PDE.)

Here comes my problem:
For general g,
if I can find a function $w\in H^1$ such that $w|_{\partial\Omega}=g$, by letting $\tilde u=u-w$, 
we have
$$
\left\{\!\!
\begin{aligned}
&L\tilde u=\tilde f,x\in\Omega\\
&\tilde u|_{\partial \Omega}=0
\end{aligned}
\right.
$$
whose solution is already known.
So how to find such a $w$?
 A: Here is an explicit construction not requiring local straightening of the boundary. 
I think it got to do the trick but I didn't check the details. Let $d(x)$ be the distance function from the point $x$ to the boundary $\partial\Omega$, denote by $K(a,t)=(4\pi t)^{-(n-1)/2}e^{-a^2/4t}$ a kernel corresponding to the fundamental solution of the heat equation with $n-1$ spatial variables. 
For $x\in \Omega$ consider the function
$$
w(x)=e^{1/4}\int_{\partial\Omega}K(|x-y|,d^2(x))g(y)ds_y.
$$
The factor $e^{1/4}$ is here to ensure that $w|_{\partial\Omega}=g$ since near the boundary in the coordinate system with axes for $x_1,\ldots,x_{n-1}$ tangent to $\partial\Omega$ at the origin we have $|x-y|^2\approx(x_1-y_1)^2+ \ldots(x_{n-1}-y_{n-1})^2+d^2(x)$. 
Using the heat kernel allows to obtain several effects. It acts as a mollifier and then $t\to+0$ it tends to the delta function. (The use of fundamental solution of the Laplace equation (or its derivatives) as a potential's kernel would't suit for this aim as  it leads to an integral equation).
But there is a problem, the potential with the heat kernel gives an anisotropic solution well described in function spaces in which the elements are to be twice as smooth with respect to spatial variables as with time variable $t$.
The idea is to use rescaling $t\to t^2$ to mend this anisotropic property, giving equal degrees for $x$ and $t$ in the heat kernel.
If $u\in H^1(\Omega)$ then $g\in H^{1/2}(\partial\Omega)$. Consider the Cauchy problem $u_t=\Delta u$ in the half-space $\mathbb R^n_+=\mathbb R^{n-1}\times(0,\infty)\ $, $u(\cdot,0)=g\in H^{1/2}(\mathbb R^{n-1})\ $. Then the solution $u$ belongs to anisotropic space $H^{3/2,3/4}(\mathbb R^n_+)$. So the derivative $u_x$ is smooth enough. The essential thing here is that is valid the following 
Lemma 
Let $g\in H^{1/2}(\mathbb R^{n-1})$ and $v(x,t)=u(x,t^2)\;$. Then $v_t\in L_2(\mathbb R^n_+)\;$.
The rest is technical. The distance function $d$ is uniformly Lipschitz in $\mathbb R^n$ so there should be no problem taking the first derivatives etc. If the boundary is smooth then $d$ is also smooth near it. So multiplying $w$ by a cutoff function of the form $\xi(d(x))$ one can get a more smooth result.
A: If $\Omega$ is a half space, you can construct w by Fourier transforms. Any proof of the inverse trace theorem in the literature will show you specifics. For general $\Omega$, you use partition of unity and local transformations which straighten out the boundary.
A: Related to Andrew's answer is the following approach.
If you know how to solve the problem with $a=1$, $b=0$, and $f=0$, you are done. Indeed, introduce $\psi\in H^1(\Omega)$ given by
$$
\Delta \psi =0 \textrm{ in }\Omega \textrm{ and } \psi=g \textrm{ on } \partial\Omega,
$$
Then your problem becomes a the variational problem for $v=u-\psi\in H^1_0(\Omega)$, namely
$$
\int_{\Omega} a \nabla v \cdot \nabla w + \int b v w = \int_{\Omega} (f-b\psi) w
- \int_{\Omega} a \nabla \psi \cdot \nabla w  
$$
for all $w\in H^1_0(\Omega)$, which is something you can plug in to your FEM solver for Dirichlet boundary conditions, since the right hand-side is nice.
Now to find $\psi$, you can either tweak you FEM solver (looking at which nodes are on the boundary etc) or solve it by an Integral Equation Method which uses the explicit Green function of the Laplacian in the free space. The advantage is that it corresponds to solving problems in $n-1$ dimensions, and if the boundary is smooth, you can expect exponential rates of convergence (not many boundary points required). Furthermore, $\psi$ is very smooth inside, and therefore can be safely evaluated on a rough grid  (again, very few points required) and interpolated in $P_1, P_2$..or whatever space you like for your FEM inside.
