Here is an example in the plane which is genuinely discontinuous up to the boundary. I try to keep it as elementary as possible. There is an if and only if criterion for continuity up to the boundary of harmonic functions, due to Wiener: this example is essentially how one shows the only if part.

Take $B_1$, and remove from it a collection of balls $\bar{B}_{r_k}(2^{-k}),$ with $k \geq 1$, centered at points on the real axis, of radii to be chosen but with $r_k \ll 2^{-k}$. Set $U_k = B_1 \setminus \bar{B}_{r_k}(2^{-k})$, and $U = (\cap_k U_k) \setminus \{0\}$ ($0$ is also removed to make $U$ open). Solve for $u_k$ satisfying $\Delta u_k = - 1$ on $U_k$ and $u_k \in H^1_0(U_k)$, and similarly $u$ with $\Delta u = -1$ on $U$ with $u \in H^1_0(U)$ (there isn't any issue with solving this using standard energy methods, even on $U$).

From the maximum principle, $u_{k+1}\leq u_k$ on $U_{k+1}$, and $u \leq u_k$ on $U$. Let us note here that $u_0$ (with no balls removed) has $u_0 \leq \frac{1}{4} < 1$ (just by directly computing $u_0 = (1 - |x|^2)/4$).

We also have that $u_k \rightarrow u$ pointwise on $U$: indeed, the $u_k$ are uniformly bounded in $H^1(U)$ by integrating by parts, $$ \int_{U_k} |\nabla u_k|^2 = \int_{U_k} u_k \leq C \|\nabla u_k\|_{L^2(U_k)}. $$ Extract a weakly convergent subsequence $u_k \rightharpoonup v \in H^1_0(B_1)$. It is straightforward to check that $v$ solves $\Delta v = -1$ on $U$, and that $v \in H^1_0(U_k)$ for every $k$. I claim (see below) that $v \in H^1_0(U)$; then from uniqueness/weak maximum principle it follows that $v = u$, that $u_k \rightarrow u$ in $L^2(U)$, and from regularity the convergence is pointwise, or even locally in $C^j$, etc.

Consider, then, $w = u_k - u_{k + 1}$. This function is harmonic on $U_{k+1}$, and vanishes on the boundary except along $\partial B_{r_{k+1}}(2^{-k})$--along which it has $w < 1$. The function $h = 1 - \frac{\log (|x - 2^{-k}|/r_k)}{\log (2/ r_k)}$ is also harmonic on $U_k$, is nonnegative on $B_1$, and is $1$ along $\partial B_{r_{k+1}}(2^{-k})$: by the maximum principle, $w \leq h$. For $x$ in the interval $[-1, 0]$, we have that $$ w \leq 1 - \frac{\log (2^{-k}/r_k)}{\log (2/r_k)} = \frac{k\log 2}{\log (2/r_k)}. $$ Now, we get to choose $r_k$ as small as we like: for example, we can choose it to be $2e^{-t 2^k}$ for a small $t$, so that $w \leq 2k t 2^{-k}$ along this entire line segment [the main point here is to choose it so that $u_k - u_{k + 1}$ is summable with small sum].

Then summing this, we see that for any $x \in [-1, 0)$, $$ |u - u_0| \leq \sum_{k = 0}^\infty |u_k - u_{k + 1}| \leq C t. $$ Choosing $t$ small enough and using that $u_0 > 1/8$ along $[-1/2, 0]$, this gives $ u > 1/16$ on this segment.

On the other hand, on each of the boundaries of balls composing $\partial U \setminus \{0\}$, we have $u = 0$. We conclude that $u$ had to be discontinuous at $0$.

To verify the claim: let $h_r = \frac{[\log|x| - \log r^2)]_+}{|\log r|}$. This function vanishes for $|x| \leq r^2$, and for $|x| \geq r$ it has $h_r \geq 1 > v$. In between, we can integrate $$ \int_{r^2 < |x| < r} |\nabla h_r|^2 = \int_{r^2}^r \frac{1}{s |\log r|^2}ds \leq \frac{C}{|\log r|}. $$ The point of this is, set $v_r = \min\{v, h_r\}$. Then $v_r \in H^1_0(U)$ (it vanishes on $B_{r^2}$, no issues near $0$ here), and $v_r = v$ outside of $B_r$. We then have $$ \int |\nabla (v_r - v)|^2 = C \int_{|x|< r} |\nabla v|^2 + |\nabla h_r|^2, $$ with the first term going to $0$ because $|\nabla v|^2$ is integrable, while the second term is controlled by $1/|\log r|\rightarrow 0$ from above. So $v_r \rightarrow v$ in $H^1(B_1)$, and $v_r \in H^1_0(U)$: it follows that $v \in H^1_0(U)$.

Here is an example in the plane which is genuinely discontinuous up to the boundary. I try to keep it as elementary as possible. There is an if and only if criterion for continuity up to the boundary of harmonic functions, due to Wiener: this example is essentially how one shows the only if part.

Take $B_1$, and remove from it a collection of balls $\bar{B}_{r_k}(2^{-k}),$ with $k \geq 1$, centered at points on the real axis, of radii to be chosen but with $r_k \ll 2^{-k}$. Set $U_k = B_1 \setminus \cup_{j = 1}^k\bar{B}_{r_j}(2^{-j})$, and $U = (\cap_k U_k) \setminus \{0\}$ ($0$ is also removed to make $U$ open). Solve for $u_k$ satisfying $\Delta u_k = - 1$ on $U_k$ and $u_k \in H^1_0(U_k)$, and similarly $u$ with $\Delta u = -1$ on $U$ with $u \in H^1_0(U)$ (there isn't any issue with solving this using standard energy methods, even on $U$).

From the maximum principle, $u_{k+1}\leq u_k$ on $U_{k+1}$, and $u \leq u_k$ on $U$. Let us note here that $u_0$ (with no balls removed) has $u_0 \leq \frac{1}{4} < 1$ (just by directly computing $u_0 = (1 - |x|^2)/4$).

We also have that $u_k \rightarrow u$ pointwise on $U$: indeed, the $u_k$ are uniformly bounded in $H^1(U)$ by integrating by parts, $$ \int_{U_k} |\nabla u_k|^2 = \int_{U_k} u_k \leq C \|\nabla u_k\|_{L^2(U_k)}. $$ Extract a weakly convergent subsequence $u_k \rightharpoonup v \in H^1_0(B_1)$. It is straightforward to check that $v$ solves $\Delta v = -1$ on $U$, and that $v \in H^1_0(U_k)$ for every $k$. I claim (see below) that $v \in H^1_0(U)$; then from uniqueness/weak maximum principle it follows that $v = u$, that $u_k \rightarrow u$ in $L^2(U)$, and from regularity the convergence is pointwise, or even locally in $C^j$, etc.

Consider, then, $w = u_k - u_{k + 1}$. This function is harmonic on $U_{k+1}$, and vanishes on the boundary except along $\partial B_{r_{k+1}}(2^{-k})$--along which it has $w < 1$. The function $h = 1 - \frac{\log (|x - 2^{-k}|/r_k)}{\log (2/ r_k)}$ is also harmonic on $U_k$, is nonnegative on $B_1$, and is $1$ along $\partial B_{r_{k+1}}(2^{-k})$: by the maximum principle, $w \leq h$. For $x$ in the interval $[-1, 0]$, we have that $$ w \leq 1 - \frac{\log (2^{-k}/r_k)}{\log (2/r_k)} = \frac{k\log 2}{\log (2/r_k)}. $$ Now, we get to choose $r_k$ as small as we like: for example, we can choose it to be $2e^{-t 2^k}$ for a small $t$, so that $w \leq 2k t 2^{-k}$ along this entire line segment [the main point here is to choose it so that $u_k - u_{k + 1}$ is summable with small sum].

Then summing this, we see that for any $x \in [-1, 0)$, $$ |u - u_0| \leq \sum_{k = 0}^\infty |u_k - u_{k + 1}| \leq C t. $$ Choosing $t$ small enough and using that $u_0 > 1/8$ along $[-1/2, 0]$, this gives $ u > 1/16$ on this segment.

On the other hand, on each of the boundaries of balls composing $\partial U \setminus \{0\}$, we have $u = 0$. We conclude that $u$ had to be discontinuous at $0$.

To verify the claim: let $h_r = \frac{[\log|x| - \log r^2)]_+}{|\log r|}$. This function vanishes for $|x| \leq r^2$, and for $|x| \geq r$ it has $h_r \geq 1 > v$. In between, we can integrate $$ \int_{r^2 < |x| < r} |\nabla h_r|^2 = \int_{r^2}^r \frac{1}{s |\log r|^2}ds \leq \frac{C}{|\log r|}. $$ The point of this is, set $v_r = \min\{v, h_r\}$. Then $v_r \in H^1_0(U)$ (it vanishes on $B_{r^2}$, no issues near $0$ here), and $v_r = v$ outside of $B_r$. We then have $$ \int |\nabla (v_r - v)|^2 = C \int_{|x|< r} |\nabla v|^2 + |\nabla h_r|^2, $$ with the first term going to $0$ because $|\nabla v|^2$ is integrable, while the second term is controlled by $1/|\log r|\rightarrow 0$ from above. So $v_r \rightarrow v$ in $H^1(B_1)$, and $v_r \in H^1_0(U)$: it follows that $v \in H^1_0(U)$.