Persistence length $L_p$ is a basic mechanical property quantifying the bending stiffness of a polymer and is defined as the relaxation of $\langle\cos(\theta)\rangle$ of bond angles of the polymer chain:
$$ \langle\cos(\theta(s))\rangle = \exp(-sl/L_p) $$
The problem of calculating $L_p$ becomes calculating $\langle\cos(\theta(s))\rangle$, if we choose the $\theta$ as the angle of adjacent bonds, i.e., $s=1$, we have:
$$ \langle\cos(\theta)\rangle = \exp(-l/L_p) $$
in simulations, the stiffness constant is given by a bending energy $U_b$, when $U_b$ is large, where excluded volume effect is negligible, we can calculate $\langle\cos(\theta)\rangle$ as
$$\langle\cos(\theta)\rangle=\frac{\int_0^\pi \cos(\theta)\sin(\theta)\exp(-\beta U_b)\mathrm{d}\theta}{\int_0^\pi\sin(\theta)\exp(-\beta U_b)\mathrm{d}\theta}$$
the $\sin(\theta)$ is a geometric weight that when 2 bonds are at an angle of $\theta$, then in 3-D space, the number of bonds is proportion to $\sin(\theta)$.
Here is an example of harmonic bending energy:
$$U_b(\theta)=-\frac{1}{2}k \theta^2$$
where $k$ is the stiffness constant in unit of $k_BT$, the $\langle\cos(\theta)\rangle$ is simply
$$\frac{e^{\frac{3}{4 k}} \left(-2 \text{erf}\left(\frac{1}{\sqrt{k}}\right)+\text{erf}\left(\frac{1+i \pi k}{\sqrt{k}}\right)+\text{erf}\left(\frac{1-i \pi k}{\sqrt{k}}\right)\right)}{2 \left(-2 \text{erf}\left(\frac{1}{2 \sqrt{k}}\right)+\text{erf}\left(\frac{1+2 i \pi k}{2 \sqrt{k}}\right)+\text{erf}\left(\frac{1-2 i \pi k}{2 \sqrt{k}}\right)\right)}$$
For large $k$ ($k\ge 6$), $\langle\cos(\theta)\rangle\simeq L(k)$, with $L(k)$ being the Langevin function.