function S = PML_function(QuadTriangle,PML)

sigma = PML.chi * PML.rho;

S = 1 + sigma/PML.Lambda^2 .* (QuadTriangle.ynodes - PML.H).^2 ...
                           .* (QuadTriangle.ynodes > PML.H);

end