Helmholtz 1D Finite Difference Approximation using Algebraic Equation

I am trying to approximate Helmholtz's wave equation using algebraic equation.

I think I have made correct algebraic equation, but it does not work; failed to approximate as figure below.

The boundary condition I used, - ux +iwu = 0 on left edge, ux + iwu = 0 on right edge.

clear
close all

%USER_PAR.m
%% Set Parameters
%%---------------------------------------------------------------
global ax bx w v h k

% define constants
nx = 200; freq = 5;
ax = 0; bx = 1;
w = 2 * pi * freq; v = 1.5; 
k = w / v; h = (bx - ax) / nx;

% define true and source equation
% -Uxx - K^2*u = S
syms x
true(x) = exp(1i*w*(x-1)) + exp(-1i*w*x) - 2;
source(x) = diff(diff(true,x));

true_u = matlabFunction(true);
source = matlabFunction(source);

clear x true

% Main 
%USER_PAR

[A, b] = algebraic_system(source, nx);

u1 = A\b; u1 = u1';
X = linspace(ax, bx, nx+1);
U = true_u(X);

E8 = norm(U(:)-u1(:),inf); E2 = norm(U(:)-u1(:),2)/sqrt(nx/4);
    fprintf('  (nx)=(%3d); (L2,L8)-error = (%.3g , %.3g)\n',nx,E2,E8);

  (nx)=(200); (L2,L8)-error = (11.8 , 9.72)

plot(U); hold on; plot(u1);

%------------------------------------------------------------------

function [A, b] = algebraic_system(source, nx)
global w h ax k

% define A, coefficient.
A = zeros(nx+1, nx+1);

A(1,1) = 2 - h^2*k^2 + 2*h*1i*w; A(1,2) = -2;
A(nx+1, nx+1) = 2 - h^2*k^2 + 2*h*1i*w; A(nx+1, nx) = -2;

for i = 2 : nx
    A(i, i-1) = -1;
    A(i, i) = 2 - h^2*k^2;
    A(i, i+1) = -1;
end

b = zeros(nx+1, 1);
for i = 1: nx+1
    b(i) = source(ax + h*(i-1));
end
b = h^2 * b;
%----------------------------------------------------------
end

[Torsten]

Use bvp4c for real and imaginary part if you have difficulties with the discretization.