Can someone please tell me why vargm() is crashing in main()?
Code:
// Q3.cpp : Defines the entry point for the console application.
//
#include <math>
#include "stdafx.h"
#include <iostream>
using namespace std;
double t, p; //azimuth and dip angles respectively
double lambda1 = 2.2;
double lambda2 = 1;//hz = 0, so this will have no influence for the 0 degrees dip case
double a11,a12,a21,a22,d11,d13,d31,d33;
double lagh,hy,hx,Y1,Y2,C0,Y,y_1,y_2,var1,var2;
double ax = 55;// range of variogram in E-W direction: fixed
double ay = 25;//range of variograms in N-S directions: fixed
double vargm();
const double pi = 3.141592653589793238462643383279502884197;
const int rows = 3;
const int cols = 3;
double lag,f_lag;
//---------------------------------------------------------
//Initializing the arrays to zero
//---------------------------------------------------------
//array that holds the dip angle, phi
double dip_array [3][3] = {0};
//array that holds the azimuth angle, theta
double azim_array [3][3] = {0};
//array that holds the anisotropy ratios
double lambda_array [3][3] = {0};
//array that holds the original lag values
double r_array [3] = {0};
//array that holds final lag values
// double fr_array [3][3] = {0};
//array that holds final lag values
double fr_array2 [3] = {0};
//---------------------------------------------------------
//Code to convert the angles to radians, acceptable in C++
//---------------------------------------------------------
//Dip Angle,p
double a(double p)
{
d11 = cos((p*pi)/180);
return d11;
}
double b(double p)
{
d13 = -1*sin((p*pi)/180);
return d13;
}
double c(double p)
{
d31 = sin((p*pi)/180);
return d31;
}
double d(double p)
{
d33 = cos((p*pi)/180);
return d33;
}
//Azimuth Angles, t
double e(double t)
{
a11 = sin((t*pi)/180);
return a11;
}
double f(double t)
{
a12 = cos((t*pi)/180);
return a12;
}
double g(double t)
{
a21 = -1*cos((t*pi)/180);
return a21;
}
double h(double t)
{
a22 = sin((t*pi)/180);
return a22;
}
void printarray (double arg[][3])
{
for (int i=0; i<rows; i++)//for each row
{
for (int j=0; j<cols; j++) //for each column value
cout << arg[i][j] << " ";
cout << "\n";
}
}
void printarray2(double arg[])
{
for (int i=0; i<rows; i++)//for each row
{
//for (int j=0; j<3; j++) //for each column value
cout << arg[i] << " ";
cout << "\n";
}
}
int main(int argc, char* argv[])
{
cout<<"ANISOTROPY MODELLING.\n";
cout<<"\nPlease enter a real number between 0 and 360 for the azimuth angle, theta, then press the 'Enter' key.\n\n";
cin>>t;
cout<<"\nPlease enter a real number between 0 and 360 for the dip angle, phi, then press the 'Enter' key.\n\n";
cin>>p;
a(p);
b(p);
c(p);
d(p);
e(t);
f(t);
g(t);
h(t);
//-----------------------------------------------------------
//Loading the matrices
//-----------------------------------------------------------
//array that holds the dip angle, phi
double dip_array [3][3] = {{d11,0,d13},{0,1,0},{d31,0,d33}};
//array that holds the azimuth angle, theta
double azim_array [3][3] = {{a11,a12,0},{a21,a22,0},{0,0,1}};
//array that holds the anisotropy ratios
double lambda_array [3][3] = {{lambda1,0,0},{0,1,0},{0,0,lambda2}};
//array that holds the original lag values
cout<<"\nPlease enter the lag distance in the E-W (x-direction) and then press the 'Enter' key.\n\n";
cin>>hx;
cout<<"\nhx = "<<hx<<endl;//" "<<"ax = " <<ax<<endl;
cout<<"\nPlease enter the lag distance in the N-S (y-direction) and then press the 'Enter' key.\n\n";
cin>>hy;
// cout<<"\nPlease enter the range in the N-S (y-direction) and then press the 'Enter' key.\n\n";
// cin>>ay;
cout<<"\nhy = "<<hy<<endl;//" "<<"ay = "<<ay<<endl;
double r_array [3] = {hx,hy,0};
//----------------------------------------------------------------
//code to check that the arrays are properly initialized in the matrices for any theta,phi
cout<<"\nThe dip array is:\n";
printarray(dip_array);
cout<<"\nThe azimuth array is:\n";
printarray(azim_array);
cout<<"\nThe anisotropy array is:\n";
printarray(lambda_array);
cout<<"\nThe original lag array is:\n";
printarray2(r_array);
//----------------------------------------------------------------
//Matrix multiplication code
for(int i=0;i<rows;i++) //for each row
{
for(int j=0;j<cols;j++) //for each column
{
for(int k=0;k<cols;k++) //for each column
{
for(int l=0;l<cols;l++) //for each column
{
fr_array2 [i] += lambda_array[i][l]*dip_array[l][k]*azim_array[k][j]*r_array[j];
}
}
}
}
//Print the final lag array
cout<<"\nThe final lag array is:\n";
printarray2(fr_array2);
//code for determining the value of the lag
for(i=0;i<rows;i++)
{
lag += pow(fr_array2[i],2);
f_lag = pow(lag,0.5);
}
//Print the final lag value
cout<<"\nThe final lag value - given the specified azimuth and dip angles - is:\n";
cout<<"lag = "<<f_lag;//Prgm crashes beyond this point!! HEEELP!!!
vargm();
//Print the final variogram value
cout<<"\n\nThe value of the variogram - given the specified azimuth and dip angles - is:\n";
cout<<"Y(h) = "<<Y<<endl;
return 0;
}
//code for spherical variograms
double vargm()
{
if (t == 0 && p == 0)
{
lagh = sqrt(pow(hx,2)+pow(hy,2));
while (lagh<ay) //lagh/ay<1
{
var1 = 1.5*sqrt(pow(lagh/ax,2));
var2 = 0.5*pow(pow(lagh/ax,2),3);
y_1 = var1 - var2;
}
while (lagh/ay>1)
{
y_1 = 1;
}
while (hx/ax<1)
{
// h = sqrt(pow(hx,2)+pow(hy,2));
var1 = 1.5*sqrt(pow(hx/ax,2));
var2 = 0.5*pow(pow(hx/ax,2),3);
y_2 = var1 - var2;
}
while (hx/ax>1)
{
y_2 = 1;
}
}
if (t != 0 || p != 0)
{
lagh = sqrt(pow(hx,2)+pow(hy,2));
while (lagh<ay)
{
hx = fr_array2[0];
hy = fr_array2[1];
var1 = 1.5*sqrt(pow((lagh/ax),2));
var2 = 0.5*pow(pow((lagh/ax),2),3);
y_1 = var1 - var2;
}
while (lagh>ay)
{
y_1 = 1;
}
while (hx/ax<1)
{
hx = fr_array2[0];
hy = fr_array2[1];
// h = sqrt(pow(hx,2)+pow(hy,2));
var1 = 1.5*sqrt(pow(hx/ax,2));
var2 = 0.5*pow(sqrt(pow(hx/ax,2)),3);
y_2 = var1 - var2;
}
while (hx/ax>1)
{
y_2 =1;
}
}
C0 = 0.125 + 0.075;
Y1 = 0.075*y_1;
Y2 = 0.15*y_2;
Y = C0 + Y1 + Y2;
return Y;
}