# set terminal png transparent nocrop enhanced font arial 8 size 420,320 
# set output 'prob.26.png'
set format x "%.1f"
set format y "%.2f"
unset key
set xzeroaxis lt 0 lw 1.000
set yzeroaxis lt 0 lw 1.000
set title "laplace (or double exponential) PDF with a = 0, b = 1" 0.000000,0.000000  font ""
set xlabel "x ->" 0.000000,0.000000  font ""
set xrange [ -5.65685 : 5.65685 ] noreverse nowriteback
set ylabel "probability density ->" 0.000000,0.000000  font ""
set yrange [ 0.00000 : 0.550000 ] noreverse nowriteback
Binv(p,q)=exp(lgamma(p+q)-lgamma(p)-lgamma(q))
arcsin(x)=invpi/sqrt(r*r-x*x)
beta(x)=Binv(p,q)*x**(p-1.0)*(1.0-x)**(q-1.0)
bin_s(x)=n!/(n-int(x))!/int(x)!*p**int(x)*(1.0-p)**(n-int(x))
bin_l(x)=exp(lgamma(n+1)-lgamma(n-int(x)+1)-lgamma(int(x)+1)+int(x)*log(p)+(n-int(x))*log(1.0-p))
binom(x)=(n<20)?bin_s(x):bin_l(x)
cauchy(x)=b/(pi*(b*b+(x-a)**2))
chi(x)=exp((0.5*df1-1.0)*log(x)-0.5*x-lgamma(0.5*df1)-df1*0.5*log2)
erlang(x)=lambda**n/(n-1)!*x**(n-1)*exp(-lambda*x)
extreme(x)=alpha*(exp(-alpha*(x-u)-exp(-alpha*(x-u))))
f(x)=Binv(0.5*df1,0.5*df2)*(df1/df2)**(0.5*df1)*x**(0.5*df1-1.0)/(1.0+df1/df2*x)**(0.5*(df1+df2))
g(x)=exp(rho*log(lambda)+(rho-1.0)*log(x)-lgamma(rho)-lambda*x)
geometric(x)=exp(log(p)+int(x)*log(1.0-p))
halfnormal(x)=sqrt2invpi/sigma*exp(-0.5*(x/sigma)**2)
hypgeo(x)=(int(x)>mm||int(x)<mm+n-nn)?0:mm!/(mm-int(x))!/int(x)!*(nn-mm)!/(n-int(x))!/(nn-mm-n+int(x))!*(nn-n)!*n!/nn!
laplace(x)=0.5/b*exp(-abs(x-a)/b)
logistic(x)=lambda*exp(-lambda*(x-a))/(1.0+exp(-lambda*(x-a)))**2
lognormal(x)=invsqrt2pi/sigma/x*exp(-0.5*((log(x)-mu)/sigma)**2)
maxwell(x)=fourinvsqrtpi*a**3*x*x*exp(-a*a*x*x)
negbin(x)=exp(lgamma(r+int(x))-lgamma(r)-lgamma(int(x)+1)+r*log(p)+int(x)*log(1.0-p))
nexp(x)=lambda*exp(-lambda*x)
normal(x)=invsqrt2pi/sigma*exp(-0.5*((x-mu)/sigma)**2)
pareto(x)=x<a?0:b/x*(a/x)**b
poisson(x)=mu**int(x)/int(x)!*exp(-mu)
rayleigh(x)=lambda*2.0*x*exp(-lambda*x*x)
sine(x)=2.0/a*sin(n*pi*x/a)**2
t(x)=Binv(0.5*df1,0.5)/sqrt(df1)*(1.0+(x*x)/df1)**(-0.5*(df1+1.0))
triangular(x)=1.0/g-abs(x-m)/(g*g)
uniform(x)=1.0/(b-a)
weibull(x)=lambda*n*x**(n-1)*exp(-lambda*x**n)
carcsin(x)=0.5+invpi*asin(x/r)
cbeta(x)=ibeta(p,q,x)
cbinom(x)=ibeta(n-x,x+1.0,1.0-p)
ccauchy(x)=0.5+invpi*atan((x-a)/b)
cchi(x)=igamma(0.5*df1,0.5*x)
cerlang(x)=1.0-exp(-lambda*x)*(1.0+lambda*x+0.5*(lambda*x)**2)
cextreme(x)=exp(-exp(-alpha*(x-u)))
cf(x)=1.0-ibeta(0.5*df2,0.5*df1,df2/(df2+df1*x))
cgamma(x)=igamma(rho,x)
cgeometric(x)=(x<1)?geometric(0):geometric(x)+cgeometric(x-1)
chalfnormal(x)=erf(x/sigma/sqrt2)
chypgeo(x)=(x<1)?hypgeo(0):hypgeo(x)+chypgeo(x-1)
claplace(x)=(x<a)?0.5*exp((x-a)/b):1.0-0.5*exp(-(x-a)/b)
clogistic(x)=1.0/(1.0+exp(-lambda*(x-a)))
clognormal(x)=cnormal(log(x))
cnormal(x)=0.5+0.5*erf((x-mu)/sigma/sqrt2)
cmaxwell(x)=igamma(1.5,a*a*x*x)
cnegbin(x)=(x<1)?negbin(0):negbin(x)+cnegbin(x-1)
cnexp(x)=1.0-exp(-lambda*x)
cpareto(x)=x<a?0:1.0-(a/x)**b
cpoisson(x)=1.0-igamma(x+1.0,mu)
crayleigh(x)=1.0-exp(-lambda*x*x)
csine(x)=x/a-sin(n*twopi*x/a)/(n*twopi)
ct(x)=(x<0.0)?0.5*ibeta(0.5*df1,0.5,df1/(df1+x*x)):1.0-0.5*ibeta(0.5*df1,0.5,df1/(df1+x*x))
ctriangular(x)=0.5+(x-m)/g-(x-m)*abs(x-m)/(2.0*g*g)
cuniform(x)=(x-a)/(b-a)
cweibull(x)=1.0-exp(-lambda*x**n)
fourinvsqrtpi = 2.25675833419103
invpi = 0.318309886183791
invsqrt2pi = 0.398942280401433
log2 = 0.693147180559945
sqrt2 = 1.4142135623731
sqrt2invpi = 0.797884560802865
twopi = 6.28318530717959
a = 0.0
alpha = 1.0
b = 1.0
df1 = 7.0
df2 = 6.0
g = 1.0
lambda = 2.0
m = 0.0
mm = 25
mu = 0.0
nn = 75
n = 10
p = 0.333333333333333
q = 2.5
r = 2.0
rho = 2.0
sigma = 1.4142135623731
u = 0.0
xmin = -5.65685424949238
xmax = 5.65685424949238
ymax = 0.55
xinc = 1
s = 0.602810274989087
plot laplace(x)