Date Class and 2D Matrix Class
Basics of Object-Oriented Programming
Problem 1:
We gave a few different versions of the DateClass.
a) Modify the setDatefunction to use a switch(m_month) statement that checks if the day is appropriate for the month given. For February assume thatthere are 29 days if the year is divisible by 4 (this is easy to check using the integer mod function(% ) using year%4 which gives the integer remainder of the year divided by 4 so if this is zero,the year is divisible by 4) and otherwise there are 28. In addition add a check to ensure that the year is in the range
of recorded history in the past (going back to ~-3,000) to some anticipated future where this
code will have passed into long obsolescence (say ~3000 to be generous).
To make testing of this function more straightforward, print out warnings rather than using assert
statements (which will stop the program). Also, have your function return a boolean value of
true if the date is successfully set and false if it does not. You should check this returned
value in the main program and print out another warning if the date was not set successfully.
b) Add a constructor to the class that accepts a month, day, and year as arguments but has a
default date of 1 Jan 2000. Your constructor should call your setDatefunction for non-default
values and use an assert statement on setDate’s return value to ensure the date was set
successfully.
c) Change the print() member function into an overload of the <<operator. Note that the
overload of the <<operator will need to be a friend function.
Problem 2:
We created an IntArray class in the notes that illustrated the concept of RAII. The
code for IntArray is posted on OWL in the folder PS3Code, IntArray.cpp.
a) Using IntArray as a template, create a new class DoubArraythat dynamically creates a
double array using the RAII paradigm. In addition to allocating the space, the constructor
should initialize the array elements to zero (You can initialize a dynamically allocated
array with uniform initialization and if you put no numbers in the brackets it is assumed
all entries are zero, e.g. array=new double[5] {}; will dynamically allocate an
array of length 5 with all entries 0).
b) Add a member function L2Norm that finds the Euclidean norm of the array assuming it
is a vector of values (i.e. if this is a vector, return its magnitude, or length). For a long
vector this can be a costly procedure so it is worthwhile to save the result as a member
variable m_normand return this value if the calculation has been previously performed
by an earlier function call. This necessitates another member variable m_normflag
that should be true if the old norm value can be reused and false if it needs to be
recalculated. These member variables should be initialized in the constructor and the
value of m_normflagshould be reset whenever the values of the array are changed.
c) Add an operator overload for the – operator and the assignment operator = that computes
the difference of two arrays as vectors (i.e. if a, b, and c are of the DoubArrayclass
then a statement like c=a-b; should be possible once these two operators have been
overloaded).
d) Further alter your code so that it uses the DoubArrayclass where the length
is dictated by the number of gridpointsN. This should be the only array type used by
your code.Write a function getNto ask the user for this value. Yourimplementation of the implicit Euler’s methods should also be turned into a functionImpEulerthat accepts as arguments N and y where N is the number of grid points andy is a DoubArrayof appropriate length passed by reference and whose memberm_array[0] is set to the initial value 1.0.
e) To test your code we will make use of the functions you created for the DoubArray
class. You need to ensure all parts of your code are tested.
Write another function to fill a DoubArraywith the exact solution (the arguments
should be the same as your ImpEulerfunction). Rather than use the getNfunction for
this test, loop over values of N=10, 20, 30, …100 and in each iteration of the loop find an
exact solution, the implicit Euler solution, the difference between the exact and Euler
solution (using the overloaded – operator), and the Euclidean norm of the difference
(using the L2Norm function). You should output to a file the value of h and the
Euclidean norm of the difference between the exact and Euler approximation and include
a plot of the norm of the difference versus h. Comment on what you find (Does the
difference appear to go to zero as h goes to zero? Is the plot linear?
Problem 3:
This problem will use the RAII concept for matrices. You are reminded that the
code for IntArray is posted on OWL in the folder PS3Code, IntArray.cpp.
a) Using the IntArray class as a template, create a class DoubMatrixthat dynamically
allocates n x m matrices. (the number of rows and number of columns should be
member variables). The matrices should be created as contiguous blocks of memory. The
copy constructor should use assertions to verify that the two matrices have the same size.
b) You have to take arbitrary sized matrices of class DoubMatrix. This should be a friend
function of the DoubMatrixclass that overloads the * operator.
c) You have to take arbitrary sized matrices of class DoubMatrixand vectors of class
DoubArrayfrom part a) of problem 2 above. The vector-matrix and matrix-vector
multiplications should be overloaded * operators that are friend functions of both
DoubMatrixand DoubArray.
As part c) now involves multiple classes, you should separate your file into a main.cpp driver
program that tests various functions, a DoubArray.cpp and DoubMatrix.cpp for the class
implementations and two header files DoubArray.h and DoubMatrix.h.
Inheritance, Virtual Functions, & Templates
Problem 1: In the last assignment, problem 2, you created the DoubArray class and added some
features in problem 2a, 2b, and 2c.
- a) Add operator overloads to the DoubArray class for the [] operator, providing both a
version for const objects and non-const objects and providing a check on the validity of
the index. The class implementation should be in a separateDoubArray.cpp file and it should havean associated header DoubArray.h.
b) we wrote several implementations of the “DList¬Series ¬Polynomial”
inheritance tree. Alter the Series and Polynomial classes as necessary to inherit from DoubArray
rather than DList You will need to addconstructors to these classes to pass up the size of the array to be initialized in theDoubList constructor.
Problem 3: Use templates to write a single class that can replace our IntArray and DoubArray
classes.
Solution
abstractodesolver.cpp
//problem 2.1
//so writing the methods
#include “abstractodesolver.h”
voidAbstractOdeSolver::SetStepSize(double h)
{
stepSize = h;
}
voidAbstractOdeSolver::SetTimeInterval(double t0, double t1)
{
initialTime = t0;
finalTime = t1;
}
voidAbstractOdeSolver::SetInitialValue(double y0)
{
initialValue = y0;
}
abstractodesolver.h
#ifndef ABSTRACTODESOLVER_H
#define ABSTRACTODESOLVER_H
classAbstractOdeSolver
{
public:
voidSetStepSize(double h);
voidSetTimeInterval(double t0, double t1);
voidSetInitialValue(double y0);
virtual double RightHandSide(double y, double t) = 0;
virtual double SolveEquation() = 0;
protected:
doublestepSize;
doubleinitialTime;
doublefinalTime;
doubleinitialValue;
};
#endif // ABSTRACTODESOLVER_H
doubarray.cpp
#include “doubarray.h”
DoubArray::DoubArray(int length) // constructor
{
assert(length > 0);
m_array = new double[length]{};
m_norm = 0.0; //because 0 init
m_normflag = true; //because 0 init
m_length = length;
//std::cout<< “…finishing constructor\n”;
}
DoubArray::DoubArray(constDoubArray&org_array) // copy constructor
{
m_normflag = org_array.m_normflag;
m_norm = org_array.m_norm;
m_length = org_array.m_length;
//delete first else memory leak !!!!
delete[] m_array;
m_array = new double[m_length];
for (int i=0; i<m_length; i++)
m_array[i]=org_array.m_array[i];
//std::cout<< “…finishing copy constructor\n”;
}
DoubArray::~DoubArray() // destructor
{
// Dynamically delete the array we allocated earlier
delete[] m_array ;
//std::cout<< “…finishing destructor\n”;
}
voidDoubArray::setValue(int index, double value)
{
reset();
m_array[index] = value;
}
voidDoubArray::reset()
{
m_normflag = false;
}
double&DoubArray::getValue(int index)
{
returnm_array[index];
}
intDoubArray::getLength()
{
returnm_length;
}
doubleDoubArray::L2Norm()
{
if(m_normflag)
returnm_norm;
m_norm = 0.0;
for(int i = 0; i <m_length; ++i)
m_norm += std::pow(m_array[i],2.0);
returnm_normflag = true, m_norm;
}
DoubArray&DoubArray::operator =(DoubArrayconst&other)
{
m_normflag = other.m_normflag;
m_norm = other.m_norm;
m_length = other.m_length;
delete[] m_array;
m_array = new double[m_length];
for (int i=0; i<m_length; i++)
m_array[i]=other.m_array[i];
return *this;
}
DoubArrayDoubArray::operator -(DoubArrayconst&other) const
{
DoubArraydoubArray(m_length);
doubArray.reset();
for (int i=0; i<m_length; i++)
doubArray.m_array[i]=m_array[i]-other.m_array[i];
returndoubArray;
}
doubleconst&DoubArray::operator [](int index)const
{
assert((index >= 0) && (m_length> index));
returnm_array[index];
}
double&DoubArray::operator [](int index)
{
assert((index >= 0) && (m_length> index));
reset();
returnm_array[index];
}
doubarray.h
//doubarray class from past assignment
#ifndef DOUBARRAY_H
#define DOUBARRAY_H
#include <iostream>
#include <cmath>
#include <cassert>
classDoubArray
{
private:
intm_length; //array length
doublem_norm; //array discrete L2Norm
boolm_normflag; //array discrete L2Norm validation flag
protected:
double *m_array; //underlying array
public:
DoubArray(int length); //constructor
DoubArray(constDoubArray&org_array); //copy constructor
virtual ~DoubArray(); //virtual destructor
virtual void reset(); //resets flags
voidsetValue(int index, double value); //set array value by index
double&getValue(int index); //get value by index
intgetLength(); //get array length
double L2Norm(); //get l2Norm
DoubArray&operator =(DoubArrayconst&other); //assignment operator overload
DoubArray operator -(DoubArrayconst&other) const; //minus operator overload
doubleconst&operator [](int index)const; //bracket op overload
double&operator [](int index); //bracket op overload
};
#endif // DOUBARRAY_H
forwardeulersolver.cpp
#include “forwardeulersolver.h”
#include <fstream>
voidForwardEulerSolver::SetRightHandSide(rightHandSide_tconst&value)
{
rightHandSide = value;
}
doubleForwardEulerSolver::RightHandSide(double y, double t)
{
returnrightHandSide(y,t);
}
voidForwardEulerSolver::SetFileName(std::string const&value)
{
fileName = value;
}
doubleForwardEulerSolver::SolveEquation()
{
//prepare file to write output
std::ofstreamfos(fileName.c_str());
if(!fos.is_open())
throwstd::runtime_error(“Can’t open file outExpEuler.data”);
//init initial value problem (IVP)
double t(initialTime);
double y(initialValue);
//write IVP to file
fos<< t << ” ” << y << “\n”;
while(t <finalTime)
{
//calculate new stepsize
double h = std::min(stepSize,finalTime-t);
//update solution
y = y + h*RightHandSide(y,t);
//update time node
t = t+h;
//write value snapshot to file
fos<< t << ” ” << y << “\n”;
}
//return final solution
return y;
}
forwardeulersolver.h
//problem 2.2
#ifndef FOWARDEULERSOLVER_H
#define FOWARDEULERSOLVER_H
#include “abstractodesolver.h”
#include <functional>
classForwardEulerSolver :
publicAbstractOdeSolver
{
public:
typedefstd::function<double(double,double)>rightHandSide_t;
doubleRightHandSide(double y, double t);
doubleSolveEquation();
//set right hand side
voidSetRightHandSide(conststd::function<double (double y, double t)>&value);
voidSetFileName(std::string const&value);
private:
rightHandSide_trightHandSide;
std::stringfileName;
};
#endif // FOWARDEULERSOLVER_H
intarray.h
#ifndef INTARRAY_H
#define INTARRAY_H
#include <iostream>
#include <cassert>
classIntArray
{
private:
int *m_array;
intm_length;
public:
IntArray(int length) // constructor
{
assert(length > 0);
m_array = new int[length];
m_length = length;
std::cout<< “…finishing constructor\n”;
}
IntArray(constIntArray&org_array) // copy constructor
{
m_length = org_array.m_length;
m_array = new int[m_length];
for (int i=0; i<m_length; i++)
m_array[i]=org_array.m_array[i];
std::cout<< “…finishing copy constructor\n”;
}
~IntArray() // destructor
{
// Dynamically delete the array we allocated earlier
delete[] m_array ;
std::cout<< “…finishing destructor\n”;
}
voidsetValue(int index, int value) { m_array[index] = value; }
int&getValue(int index) { return m_array[index]; }
intgetLength() { return m_length; }
};
#endif // INTARRAY_H
laplacesolver.h
//able to calculate potential (EvalSum) given radius r and angle theta
#ifndef LAPLACESOLVER_H
#define LAPLACESOLVER_H
#include “rthetalegendre.h”
classLaplaceSolver
{
public:
typedef double(*V_t)(double);
private:
RThetaLegendre *myRTheta; //polar legendre class made in previous task
double *angles; //array containing theta’s
doubledTheta; //stepsize in theta
int M; //angles size
V_t V; //pased user defined function
public:
//constructor
LaplaceSolver(int N,
int M,
V_t V)
{
//allocmemmyRTheta
myRTheta = new RThetaLegendre(N);
//discrete version of coefficients is not correct in task 1e)
//integral should be approxed from 0 to pi with M points
//so this gives stepsize h[i] = i*pi/(M-1) so h[0] = 0 and h[M-1] = pi
//1) make size M array
angles = new double[M];
//stepsizedTheta is pi/(size(angles)-1) !!!
dTheta = M_PI/(double(M-1));
//populate angles correct formula theta = i*pi/(M-1) = i*dTheta
for(int i = 0; i < M; ++i)
angles[i] = double(i)*dTheta;
this->M = M;
this->V = V;
}
~LaplaceSolver()
{
deletemyRTheta;
delete[] angles;
}
//evaluate sum
doubleEvalSum(double r, double theta)
{
returnmyRTheta->EvalSum(r,theta);
}
//evaluates the coefficients
voidEvalCoefficients()
{
//bluff because r^(-n-1) = 1.0 if r = 1.0
myRTheta->SetR(1.0);
for(int i = 0; i < M; ++i) //start with integral loop in order to keep x valid in polynomial calc
{
double theta = angles[i];
doublesinTheta = std::sin(theta);
doublecosTheta = std::cos(theta);
doubleVTheta = V(theta);
doubletmp = VTheta*sinTheta*dTheta;
for(int n = 0; n <myRTheta->getLength(); ++n)
(*myRTheta)[n] += tmp*(n+1.0/2.0)*myRTheta->EvalPoly(n,cosTheta);
}
}
};
#endif // LAPLACESOLVER_H
legendre.h
#ifndef LEGENDRE_H
#define LEGENDRE_H
#include “series.h”
class Legendre :
public Series
{
public:
doublem_lft;
doublem_rgt;
doublem_h;
Legendre(int N,
doublem_lft = -1.0,
doublem_rgt = +1.0) :
Series(N)
{
this->m_lft = m_lft;
this->m_rgt = m_rgt;
m_h = (this->m_rgt – this->m_lft)/double(N-1);
}
virtual ~Legendre()
{
}
//evaluate polynomial (warning x must stay same until reset triggered with n = 0)
virtual double EvalPoly(int n, double x)
{
if (!((x <= m_rgt) && (x >=m_lft)))
std::cout<< “Warning: Outside expected range for this polynomial. \n”;
static double Pnm1, Pn, myX;
if(n == 0)
{
myX = x;
returnPn = 1.0;
}
else if(n > 0)
{
if(myX != x)
{
std::cout<< “warning x must stay same until reset triggered with n = 0”;
exit(1);
}
//substite n -> n-1 from sheet formula
Pnm1 = ((2.0*n-1.0)*x*Pn – (n-1.0)*Pnm1)/double(n);
//swap pnm1, pn
std::swap(Pnm1,Pn);
returnPn;
}
//not allowed to reach if so return nan
return 1.0/0.0;
}
doubleEvalSum(double x)
{
doubletmp(0.0);
for(int i = 0; i <getLength(); ++i)
tmp+=m_array[i]*EvalPoly(i,x);
returntmp;
}
};
#endif // LEGENDRE_H
main.cpp
#include <iostream>
#include <fstream>
#include <cmath>
#include “intarray.h”
#include “doubarray.h”
#include “marray.h”
#include “polynomial.h”
#include “legendre.h”
#include “rthetalegendre.h”
#include “laplacesolver.h”
#include <sstream>
#include “forwardeulersolver.h”
#include “rungekuttasolver.h”
//problem 1a)
voiddoubArrayBracketOperatorDemo()
{
//test by simply check if c = a-b
DoubArraya(1), b(1), c(1);
a[0] = 5.0;
b[0] = 3.0;
c[0] = a[0]-b[0];
assert(c[0] == 2.0);
std::cout<< “Problem 1a test passed” <<std::endl;
}
//problem 1b)
voidinheritDoubArrayDemo()
{
//analogous test as for problem 1a)
{
Polynomial a(1), b(1), c(1);
a[0] = 5.0;
b[0] = 3.0;
c[0] = a[0]-b[0];
assert(c[0] == 2.0);
}
//analogous test as for problem 1a)
{
Series a(1), b(1), c(1);
a[0] = 5.0;
b[0] = 3.0;
c[0] = a[0]-b[0];
assert(c[0] == 2.0);
}
std::cout<< “Problem 1b test passed” <<std::endl;
}
//problem 1c)
voidlegendrePolyDemo()
{
//open for output
std::ofstreamfos(“legendre5.dat”);
//make Legendre series with 30 polynomials
Legendre legendre(30);
//evaluate polynomials
for(double x = legendre.m_lft; x <legendre.m_rgt; x+= legendre.m_h)
{
fos<< x << ” “;
for(int i = 0; i < 5; ++i)
{
legendre[i] = 1.0;
fos<<legendre.EvalSum(x) << ” “;
legendre[i] = 0.0;
}
fos<< “\n”;
}
std::cout<< “Problem 1c legendre5.dat written (plot by e.g. gnuplot plot.gp)” <<std::endl;
}
//problem 1d)
voidrThetaLegendreDemo()
//an = 1/n not possible due to singularity @ n = 0
//so i chose an = 1/n, and a0 = 1
{
//open for output
std::ofstreamfos(“rThetaLegendre4.dat”);
//set an = 1/(n+1)
RThetaLegendrerThetaLegendre(100);
rThetaLegendre[0] = 1.0;
for(int i = 1; i <rThetaLegendre.getLength(); ++i)
rThetaLegendre[i] = 1.0/double(i);
//evaluate series
for(double theta = 0.0; theta < M_PI; theta+= M_PI/99.0)
{
fos<< theta << ” “;
fos<<rThetaLegendre.EvalSum(1,theta) << ” “;
fos<<rThetaLegendre.EvalSum(2,theta) << ” “;
fos<<rThetaLegendre.EvalSum(4,theta) << ” “;
fos<<rThetaLegendre.EvalSum(8,theta) << ” “;
fos<< “\n”;
}
std::cout<< “Problem 1d rThetaLegendre4.dat written (plot by e.g. gnuplot plot.gp)” <<std::endl;
}
//problem 1e)
voidlaplaceSolverDemo()
{
//open for output
std::ofstreamfos(“laplaceSolver4.dat”);
LaplaceSolverlaplaceSolver(/*evenly spaced points [0,pi]*/100,
/*numPoints for integral approximation*/1e4,
/*V(theta)=*/[](double theta){ return std::sin(2.0*theta);});
//set coefficients
laplaceSolver.EvalCoefficients();
//evaluate series
for(double theta = 0.0; theta < M_PI; theta+= M_PI/99.0)
{
fos<< theta << ” “;
fos<<laplaceSolver.EvalSum(1,theta) << ” “;
fos<<laplaceSolver.EvalSum(2,theta) << ” “;
fos<<laplaceSolver.EvalSum(4,theta) << ” “;
fos<<laplaceSolver.EvalSum(8,theta) << ” “;
fos<< “\n”;
}
std::cout<< “Problem 1e laplaceSolver4.dat written (plot by e.g. gnuplot plot.gp)” <<std::endl;
}
//problem 2)
void task3TestSolver()
{
//lambda the RightHandSide given
autoRightHandSide = [](double /*y*/, double t) { return 1.0 + t;};
//lambda the exact solution given
autoExactSolution = [](double t){ return 0.5*(t*t+2*t+4.0);};
//init explicit euler
ForwardEulerSolverforwardEulerSolver;
RungeKuttaSolverrungeKuttaSolver;
//initial values for IVP
doubleinitialTime(0.0);
doublefinalTime(1.0);
doubleinitialValue(2.0);
doublestepSize(0.1);
//init solver
forwardEulerSolver.SetInitialValue(initialValue);
forwardEulerSolver.SetRightHandSide(RightHandSide);
forwardEulerSolver.SetTimeInterval(initialTime,finalTime);
forwardEulerSolver.SetStepSize(stepSize);
rungeKuttaSolver.SetInitialValue(initialValue);
rungeKuttaSolver.SetRightHandSide(RightHandSide);
rungeKuttaSolver.SetTimeInterval(initialTime,finalTime);
rungeKuttaSolver.SetStepSize(stepSize);
//1) check solution, by checking relative error
// ||yApprox – yExact||/||yExact|| <rtol
{
//set output file name
forwardEulerSolver.SetFileName(“ForwardEulerSolver_h01Run.dat”);
rungeKuttaSolver.SetFileName(“RungeKuttaSolver_h01Run.dat”);
double yApprox1 = forwardEulerSolver.SolveEquation();
std::cout<< “ForwardEulerSolver_h01Run.dat written” <<std::endl;
double yApprox2 = rungeKuttaSolver.SolveEquation();
std::cout<< “RungeKuttaSolver_h01Run.dat written” <<std::endl;
double relErr1 = std::abs(yApprox1 – ExactSolution(1.0))/std::abs(ExactSolution(1.0));
double relErr2 = std::abs(yApprox2 – ExactSolution(1.0))/std::abs(ExactSolution(1.0));
std::cout<< “ForwardEulerSolver with stepsize 0.001 has relative error = ” << relErr1 <<std::endl;
std::cout<< “RungeKuttaSolver with stepsize 0.001 has relative error = ” << relErr2 <<std::endl;
}
//2) how the choice of the size affects the accuracy ?
// error goes like E ~ h ^ p, with E:error, h:stepSize and p order of method
// so E1/E2 = (h1/h2)^p ~> p = log(E1/E2)/log(h1/h2)
// now expecting p = 1 for euler and p = 4 for rk4
// but rightHandSide only scales with t so we have exact solution with all integrators where p > 1
// therefore for rk4 we have independant of stepSize the exact solution in this case (apart from possible instabilities)
{
//estimate order for euler method
double h1 = 0.5*(finalTime-initialTime);
forwardEulerSolver.SetStepSize(h1);
forwardEulerSolver.SetFileName(“ForwardEulerSolver_h05Run.dat”);
double yApprox1 = forwardEulerSolver.SolveEquation();
std::cout<< “ForwardEulerSolver_h05Run.dat written” <<std::endl;
double h2 = (finalTime-initialTime);
forwardEulerSolver.SetStepSize(h2);
forwardEulerSolver.SetFileName(“ForwardEulerSolver_h10Run.dat”);
double yApprox2 = forwardEulerSolver.SolveEquation();
std::cout<< “ForwardEulerSolver_h10Run.dat written” <<std::endl;
double p = std::log(std::abs(yApprox1-ExactSolution(finalTime))/std::abs(yApprox2-ExactSolution(finalTime)))/std::log(h1/h2);
std::cout<< “Estimated ForwardEulerSolver order is “<< p <<std::endl;
}
}
//problem 3)
voidmArrayDemo()
{
//test by comparing sum of Array<int> Array<double>
int n = 50;
Array<int>iArray(n);
Array<double>dArray(n);
//initi
for(int i = 0; i < n; ++i)
{
iArray[i] = 4;
dArray[i] = 4.5;
}
//calculate sum
intiSum(0);
doubledSum(0);
for(int i = 0; i < n; ++i)
{
iSum += iArray[i];
dSum += dArray[i];
}
//check if correct
assert(iSum == 50*4);
assert(std::abs(dSum – 4.5*50.0) < 1.e-6);
std::cout<< “Problem 3 test passed” <<std::endl;
}
int main(int , char *[])
{
std::cout<< “The beginning” <<std::endl;
std::cout<< “\nProblem1” <<std::endl;
doubArrayBracketOperatorDemo();
inheritDoubArrayDemo();
legendrePolyDemo();
rThetaLegendreDemo();
laplaceSolverDemo();
//problem2)
std::cout<< “\nProblem2” <<std::endl;
task3TestSolver();
std::cout<< “\nProblem3” <<std::endl;
mArrayDemo();
std::cout<< “\nThe end” <<std::endl;
return 0;
}
marray.h
//Problem 3)
//basically replaced all double inside DoubArray with template argument _T
//keep in mind nobody would use these arrays in real life…
//nowadays expression templates are used in state of the art implementations
#ifndef MARRAY_H
#define MARRAY_H
template<typename _T>
class Array
{
private:
intm_length;
_T m_norm;
boolm_normflag;
protected:
_T *m_array;
public:
Array(int length) // constructor
{
assert(length > 0);
m_array = new _T[length]{};
m_norm = 0.0; //because 0 init
m_normflag = true; //because 0 init
m_length = length;
//std::cout<< “…finishing constructor\n”;
}
Array(const Array &org_array) // copy constructor
{
m_normflag = org_array.m_normflag;
m_norm = org_array.m_norm;
m_length = org_array.m_length;
//delete first else memory leak !!!!
delete[] m_array;
m_array = new _T[m_length];
for (int i=0; i<m_length; i++)
m_array[i]=org_array.m_array[i];
//std::cout<< “…finishing copy constructor\n”;
}
~Array() // destructor
{
// Dynamically delete the array we allocated earlier
delete[] m_array ;
//std::cout<< “…finishing destructor\n”;
}
virtual void reset()
{
m_normflag = false;
}
voidsetValue(int index, _T value)
{
reset();
m_array[index] = value;
}
_T&getValue(int index)
{
returnm_array[index];
}
intgetLength()
{
returnm_length;
}
_T L2Norm()
{
if(m_normflag)
returnm_norm;
m_norm = 0.0;
for(int i = 0; i <m_length; ++i)
m_norm += std::pow(m_array[i],2.0);
returnm_normflag = true, m_norm;
}
Array &operator =(Array const&other)
{
m_normflag = other.m_normflag;
m_norm = other.m_norm;
m_length = other.m_length;
delete[] m_array;
m_array = new _T[m_length];
for (int i=0; i<m_length; i++)
m_array[i]=other.m_array[i];
return *this;
}
Array operator -(Array const&other)
{
Array<_T>mArray(m_length);
mArray.reset();
for (int i=0; i<m_length; i++)
mArray.m_array[i]=m_array[i]-other.m_array[i];
returnmArray;
}
_T const&operator [](int index)const
{
assert((index >= 0) && (m_length> index));
returnm_array[index];
}
_T &operator [](int index)
{
assert((index >= 0) && (m_length> index));
reset();
returnm_array[index];
}
};
#endif // MARRAY_H
polynomial.h
#ifndef POLYNOMIAL_H
#define POLYNOMIAL_H
#include <iostream>
#include <cassert>
#include <cmath>
#include “series.h”
class Polynomial :
public Series
{
private :
doublem_lft=-1.;
doublem_rgt=1.;
public:
Polynomial(double leftendpoint=-1.0,
doublerightendpoint=1.0,
int N=10) : Series(N)
{
m_lft=leftendpoint;
m_rgt=rightendpoint;
//std::cout<< “…constructing a Polynomial” <<std::endl;
}
doubleEvalPoly(double x)
{
if (!((x <m_rgt) && (x>=m_lft)))
std::cout<< x << “Warning: Outside range for this polynomial. \n”;
doublepx=0., xn=1.0;
for (int n = 0; n <getLength(); n++)
{
px += m_array[n]*xn;
xn *= x;
}
returnpx;
}
virtual ~Polynomial()
{
}
bool decreasing(const double x)
{
intmaxN=getLength();
doublexn=1.0;
for (int i = 1; i <maxN; i++)
if (fabs(m_array[i-1]*xn) <fabs(m_array[i]*(xn*x)))
{
xn=xn*x;
return false;
}
return true;
}
};
#endif // POLYNOMIAL_H
rthetalegendre.h
#ifndef RTHETALEGENDRE_H
#define RTHETALEGENDRE_H
#include “legendre.h”
classRThetaLegendre :
//legendre base
public Legendre
{
private:
//use hint and make r member
double r;
public:
//use base constructor
using Legendre::Legendre;
//make virtual destructor
virtual ~RThetaLegendre() { }
//Evaluate polynoms
doubleEvalPoly(int n, double cosTheta)
//polynoms are pow(r,-(n+1))*Pn(cos(theta))
{
returnstd::pow(r,-(n+1.0))*Legendre::EvalPoly(n,cosTheta);
}
//set r
voidSetR(double value)
{
this->r = value;
}
//evaluate sum
doubleEvalSum(double r, double theta)
{
//set r
SetR(r);
//use base’s evalsum
return Legendre::EvalSum(std::cos(theta));
}
};
#endif // RTHETALEGENDRE_H
rungekuttasolver.h
//problem 2
#ifndef RUNGEKUTTASOLVER_H
#define RUNGEKUTTASOLVER_H
#include “abstractodesolver.h”
#include <functional>
classRungeKuttaSolver :
publicAbstractOdeSolver
{
public:
typedefstd::function<double(double,double)>rightHandSide_t;
doubleRightHandSide(double y, double t);
doubleSolveEquation();
voidSetRightHandSide(conststd::function<double (double y, double t)>&value);
voidSetFileName(std::string const&value);
private:
rightHandSide_trightHandSide;
std::stringfileName;
};
#endif // RUNGEKUTTASOLVER_H
series.h
#ifndef SERIES_H
#define SERIES_H
#include <iostream>
#include <cassert>
#include <cmath>
#include “doubarray.h”
class Series :
publicDoubArray
{
private :
doublem_sum;
boolm_sumflag;
public:
Series(int length) :
DoubArray(length)
{
m_sum = 0.0;
reset();
}
virtual ~Series()
{
}
virtual void reset()
{
DoubArray::reset(); //first reset base
m_sumflag = false;
}
bool decreasing()
{
intmaxN=getLength();
for (int i = 1; i <maxN; i++)
if (fabs(m_array[i]) >= fabs(m_array[i-1]))
return false;
return true;
}
doublegetSum()
{
if (m_sumflag)
returnm_sum;
m_sum=0.0;
intmaxN=getLength();
for (int i =0; i <maxN; i++)
m_sum += m_array[i];
m_sumflag = true;
returnm_sum;
}
};
#endif // SERIES_H