code: public static RandomVariable bachelierOptionValue(
  			RandomVariable forward,
  			RandomVariable volatility,
  			double optionMaturity,
  			double optionStrike,
  			RandomVariable payoffUnit)
  	{
  		if(optionMaturity < 0) {
  			return forward.mult(0.0);
  		}
  		else
  		{
  			RandomVariable integratedVolatility = volatility.mult(Math.sqrt(optionMaturity));
  			RandomVariable dPlus	= forward.sub(optionStrike).div(integratedVolatility);
  
  			RandomVariable valueAnalytic = dPlus.apply(NormalDistribution::cumulativeDistribution).mult(forward.sub(optionStrike))
  					.add(dPlus.apply(NormalDistribution::density).mult(integratedVolatility)).mult(payoffUnit);
  
  			return valueAnalytic;
  		}
  	}
docs: Calculates the option value of a call, i.e., the payoff max(S(T)-K,0) P, where S follows a
  normal process with constant volatility, i.e., a Bachelier model
  \[
  \mathrm{d} S(t) = r S(t) \mathrm{d} t + \sigma \mathrm{d}W(t)
  \]
  
  @param forward The forward of the underlying \( F = S(T) \exp(r T) \).
  @param volatility The Bachelier volatility \( \sigma \).
  @param optionMaturity The option maturity T.
  @param optionStrike The option strike.
  @param payoffUnit The payoff unit (e.g., the discount factor)
  @return Returns the value of a European call option under the Bachelier model.