clear; clc;

S0 = 100;
K = 110;
r = 0.01;
v = 0.25;
T = 1;
n = 1e4;

S = S0 * exp((r - v^2 /2) * T + v * sqrt(T) * randn(1, n));
S = max(S - K, 0);
c = sum(S)/ n * exp(-r * T)

d1 = (log(S0/K) + (r + v^2))/ v /sqrt(T);
d2 = d1 - v * sqrt(T);
c_ = normcdf(d1, 0, 1) * S0 - normcdf(d2, 0, 1) * K * exp(-r*T)