p2p.wrox.com Forums

Need to download code?

View our list of code downloads.


  Return to Index  

c_plus_plus_programming thread: program help


Message #1 by "Ashley Metz" <metzal@j...> on Tue, 19 Feb 2002 15:48:49
This C++ program is supposed to price an up and out barrier call option, 

using antithetic acceleration.  It runs, but gives wrong answers.  If 

anyone can figure out what I'm doing wrong, I would greatly appreciate 

it.  Thanks!!! 



#include<iostream.h>

#include<math.h>



float ran0(int &idum);

float gasdev(int &idum);



// NOTE: This program uses antithetic acceleration



int main()

{



    float S=53;

	float X=60;

	float r=.05;

	float T=.4;

	float div=0.025;

    float upbar=75;

	float normhold;

	

	float alpha=.09;

	float beta=3.29;

	float sigma=.26;

	float rho=-.79;



	int seed=12345, ctr1, ctr2;

	

	const int M=50;

    const int N1=100;



//NOTE: N1 MUST be an even number for this program to function properly



		

	float x1, x2, E1, E2, var1, var2;

	

	float st1,st2, deltt=T/static_cast<float>(M), payoff1,payoff2;



	float callval=0;

	float gy[N1], s_hat=0, lowbound=0, upbound=0;



    for (ctr1=0;ctr1<N1/2;ctr1++)

	{ 

		st1=S;

		st2=S;

		//cout << ctr1<<"\n";

		ctr2=0;

		while(ctr2<M)

		//for (ctr2=0; ctr2<M; ctr2++)

		{

	     

			x1=gasdev(seed);

			x2=gasdev(seed);



			E1=x1;

			E2=rho*x1+(sqrt(1-pow(rho,2)))*x2;





			normhold=gasdev(seed);



			//st=st*exp( (r-div-.5*(pow(vol,2)))*deltt  +  

vol*sqrt(deltt)*gasdev(seed)  );

            

			st1= st1+ (r-div)*deltt+sqrt(var1)*st1*sqrt(deltt)

*E1;

			st2= st2+ (r-div)*deltt+sqrt(var2)*st2*sqrt(deltt)*

(-E1);

			

			var1= var1+ (alpha-beta*var1)*deltt+sigma*sqrt

(var1*deltt)*E2;

			var2= var2+ (alpha-beta*var2)*deltt+sigma*sqrt

(var2*deltt)*(-1*E2);

			

		



			if(st1>upbar)

			{

			    st1=0;

			}



			if(st2>upbar)

			{

				st2=0;

			}



			 if ((st1 ==0) && (st2==0))

			 {	

				 ctr2=M;

			 }





			ctr2++;



		}

			

		payoff1=0;

		payoff2=0;



		if (st1-X > 0)

		{

			payoff1 = (st1-X)*exp(-r*T);

		}



		if (st2-X>0)

		{

			payoff2 = (st2-X)*exp(-r*T);

		}



        gy[ctr1]=payoff1;

		gy[ctr1+(N1/2)]=payoff2;

		callval=callval+(1/static_cast<float>(N1))*payoff1+

(1/static_cast<float>(N1))*payoff2;

		//cout << ctr1 << "  "<<ctr1+(N1/2) << "\n";





	}



	for (ctr1=0; ctr1<N1; ctr1++)

	{

		s_hat=s_hat+(1/(static_cast<float>(N1)))*pow((gy[ctr1]-

callval),2);



	}















	lowbound= callval - (1.96*sqrt(s_hat))/ sqrt(static_cast<float>

(N1));

	upbound = callval + (1.96*sqrt(s_hat))/ sqrt(static_cast<float>

(N1));



	cout <<"Simulation option value:  "<< callval << "\n";



	cout <<"Simulation estimate standard deviation: "<<sqrt(s_hat)/sqrt

(static_cast<float>(N1));

    cout << "\n";



	cout << "The upper 95% condifence interval is: "<< lowbound 

<< "\n";





	cout << "The lower 95% condifence interval is: " << upbound 

<< "\n";









	return 0;

}





float ran0(int &idum)

{

    int k;



	const int IA = 16807, IM = 2147483647, IQ = 127773;

	const int IR = 2836;

	const float AM = 1/static_cast<float>(IM);



// Notice that, although idum is created as a reference with &, the syntax 

is

//  such one does not put the "&" in from of the variable during 

calculations.



    k=idum/IQ;

	idum = IA*(idum-k*IQ)-IR*k;



	if (idum < k )

	{

		idum=idum+IM;

	}



	return AM*idum;

}



float gasdev(int &idum)

{

// This function uses the Box_Mueller method for generating standard-

normal deviates.



	float fac,v1=0,v2=0, rsq=0;





		while ((rsq>=1) || (rsq ==0) )

		{

			v1=2*ran0(idum)-1;

			v2=2*ran0(idum)-1;

			rsq=pow(v1,2)+pow(v2,2);



		}





		fac = sqrt(-2*log(rsq)/rsq);





//	    cout << "fac =\t"<<fac << "\nv2= " << v2 << "\n";





	return v2*fac;

}  
Message #2 by "John Elkins" <john@v...> on Tue, 19 Feb 2002 11:29:14 -0500
Ashley:



Send me the Fortran program you translated this from and I'll try to help

you out.



j



John Elkins

Web and Database Technologies.  Storage Systems

Vermont Database Corporation

400 Upper Hollow Hill Road

Stowe VT  05672-4510 USA

802-249-0914;  xxx-xxx-xxxx  (FAX);  xxx-xxx-xxxx  (residence)

john@v... <mailto:john@v...>

www.vermontdatabase.com <http://www.vermontdatabase.com>





> -----Original Message-----

> From: Ashley Metz [mailto:metzal@j...]

> Sent: Tuesday, February 19, 2002 3:49 PM

> To: C++_Programming

> Subject: [c_plus_plus_programming] program help

>

>

> This C++ program is supposed to price an up and out barrier call option,

> using antithetic acceleration.  It runs, but gives wrong answers.  If

> anyone can figure out what I'm doing wrong, I would greatly appreciate

> it.  Thanks!!!

>

> #include<iostream.h>

> #include<math.h>

>

> float ran0(int &idum);

> float gasdev(int &idum);

>

> // NOTE: This program uses antithetic acceleration

>

> int main()

> {

>

>     float S=53;

> 	float X=60;

> 	float r=.05;

> 	float T=.4;

> 	float div=0.025;

>     float upbar=75;

> 	float normhold;

>

> 	float alpha=.09;

> 	float beta=3.29;

> 	float sigma=.26;

> 	float rho=-.79;

>

> 	int seed=12345, ctr1, ctr2;

>

> 	const int M=50;

>     const int N1=100;

>

> //NOTE: N1 MUST be an even number for this program to function properly

>

>

> 	float x1, x2, E1, E2, var1, var2;

>

> 	float st1,st2, deltt=T/static_cast<float>(M), payoff1,payoff2;

>

> 	float callval=0;

> 	float gy[N1], s_hat=0, lowbound=0, upbound=0;

>

>     for (ctr1=0;ctr1<N1/2;ctr1++)

> 	{

> 		st1=S;

> 		st2=S;

> 		//cout << ctr1<<"\n";

> 		ctr2=0;

> 		while(ctr2<M)

> 		//for (ctr2=0; ctr2<M; ctr2++)

> 		{

>

> 			x1=gasdev(seed);

> 			x2=gasdev(seed);

>

> 			E1=x1;

> 			E2=rho*x1+(sqrt(1-pow(rho,2)))*x2;

>

>

> 			normhold=gasdev(seed);

>

> 			//st=st*exp( (r-div-.5*(pow(vol,2)))*deltt  +

> vol*sqrt(deltt)*gasdev(seed)  );

>

> 			st1= st1+ (r-div)*deltt+sqrt(var1)*st1*sqrt(deltt)

> *E1;

> 			st2= st2+ (r-div)*deltt+sqrt(var2)*st2*sqrt(deltt)*

> (-E1);

>

> 			var1= var1+ (alpha-beta*var1)*deltt+sigma*sqrt

> (var1*deltt)*E2;

> 			var2= var2+ (alpha-beta*var2)*deltt+sigma*sqrt

> (var2*deltt)*(-1*E2);

>

>

>

> 			if(st1>upbar)

> 			{

> 			    st1=0;

> 			}

>

> 			if(st2>upbar)

> 			{

> 				st2=0;

> 			}

>

> 			 if ((st1 ==0) && (st2==0))

> 			 {

> 				 ctr2=M;

> 			 }

>

>

> 			ctr2++;

>

> 		}

>

> 		payoff1=0;

> 		payoff2=0;

>

> 		if (st1-X > 0)

> 		{

> 			payoff1 = (st1-X)*exp(-r*T);

> 		}

>

> 		if (st2-X>0)

> 		{

> 			payoff2 = (st2-X)*exp(-r*T);

> 		}

>

>         gy[ctr1]=payoff1;

> 		gy[ctr1+(N1/2)]=payoff2;

> 		callval=callval+(1/static_cast<float>(N1))*payoff1+

> (1/static_cast<float>(N1))*payoff2;

> 		//cout << ctr1 << "  "<<ctr1+(N1/2) << "\n";

>

>

> 	}

>

> 	for (ctr1=0; ctr1<N1; ctr1++)

> 	{

> 		s_hat=s_hat+(1/(static_cast<float>(N1)))*pow((gy[ctr1]-

> callval),2);

>

> 	}

>

>

>

>

>

>

>

> 	lowbound= callval - (1.96*sqrt(s_hat))/ sqrt(static_cast<float>

> (N1));

> 	upbound = callval + (1.96*sqrt(s_hat))/ sqrt(static_cast<float>

> (N1));

>

> 	cout <<"Simulation option value:  "<< callval << "\n";

>

> 	cout <<"Simulation estimate standard deviation: "<<sqrt(s_hat)/sqrt

> (static_cast<float>(N1));

>     cout << "\n";

>

> 	cout << "The upper 95% condifence interval is: "<< lowbound

> << "\n";

>

>

> 	cout << "The lower 95% condifence interval is: " << upbound

> << "\n";

>

>

>

>

> 	return 0;

> }

>

>

> float ran0(int &idum)

> {

>     int k;

>

> 	const int IA = 16807, IM = 2147483647, IQ = 127773;

> 	const int IR = 2836;

> 	const float AM = 1/static_cast<float>(IM);

>

> // Notice that, although idum is created as a reference with &,

> the syntax

> is

> //  such one does not put the "&" in from of the variable during

> calculations.

>

>     k=idum/IQ;

> 	idum = IA*(idum-k*IQ)-IR*k;

>

> 	if (idum < k )

> 	{

> 		idum=idum+IM;

> 	}

>

> 	return AM*idum;

> }

>

> float gasdev(int &idum)

> {

> // This function uses the Box_Mueller method for generating standard-

> normal deviates.

>

> 	float fac,v1=0,v2=0, rsq=0;

>

>

> 		while ((rsq>=1) || (rsq ==0) )

> 		{

> 			v1=2*ran0(idum)-1;

> 			v2=2*ran0(idum)-1;

> 			rsq=pow(v1,2)+pow(v2,2);

>

> 		}

>

>

> 		fac = sqrt(-2*log(rsq)/rsq);

>

>

> //	    cout << "fac =\t"<<fac << "\nv2= " << v2 << "\n";

>

>

> 	return v2*fac;

> }




> $subst('Email.Unsub').




  Return to Index