- #1
- 53
- 0
I have created a matrix with a class called Lattice. The lattice is filled with objects of type 'Dipole' which is created with another class. The problem I am having is with boundary conditions when I look for a neighbour. e.g If i pick a dipole on the top row, I want its 'above' neighbour to be the site on the bottom row etc...
Below is the .cpp and .h file which contains the functions im having problems with. Note that the value passed to both functions can be negative. I thought the %xSize should take care of this but still running into problems. I believe the dipole functions are working correctly.
EDIT:: When it crashes it is always at an edge of the matrix. Usually the top row.
If my dipole is top row, it crashes when I find the 'above' neighbour.
When returning properties of this neighbour, "type" and "direction," defined in the dipole class to be 0,1 or 0,1,2,3 it returns large +ve or -ve values for the neighbour
PLEASE help as I cant see any problems :S
Below is the .cpp and .h file which contains the functions im having problems with. Note that the value passed to both functions can be negative. I thought the %xSize should take care of this but still running into problems. I believe the dipole functions are working correctly.
Code:
#ifndef lattice
#define lattice
#include "dipole.h"
class Lattice{
private:
Dipole* dipoleArray;
double E; // pointing in +ve y axis V/m
public:
const int xSize;
const int ySize;
Lattice(const int, const int);
~Lattice();
Dipole* GetDipole(int x, int y); // return a dipole
double calcLocalEnergy( int x, int y);
};
#endif
///////////////////////////////////////////////////////////
#include "dipole.h"
#include "lattice.h"
#include <cmath>
#include <cstdlib>
#include<iostream>
Lattice::Lattice(const int xSize, const int ySize): xSize(xSize), ySize(ySize){
dipoleArray = new Dipole[xSize*ySize]; //store 1D array list
for(int i=0; i<xSize; i++){
for(int j=0; j<ySize; j++){
//loop through the entire array and create dipoles
dipoleArray[ i+j*xSize ] = Dipole(); // i + j*xSize to mimic 2D lattice form
}
}
}
Dipole* Lattice::GetDipole(int x, int y){
return &dipoleArray[ (x%xSize) + (y%ySize)*xSize ]; // modulo is for loop behaviour
}
double Lattice::calcLocalEnergy( int x, int y){
double k = 1/(4*3.14158*8.85E-12); double p1=0; double p2 =0; double r = 1E-8; double energy = 0;
double E =1e4; double pZero = 1e-27; double pOne = 1e-29;
// std::cout<<"0"<<std::endl;
Dipole* central = GetDipole(x,y); // define central dipole in the lattice at position [x,y] given by the argument
// take modulo because x,y can potentially = -1 OR xSize for neighbour dipoles
// modulo ensures neighbours are found correctly
// for 0 <= x,y < xSize modulo has no effect
Dipole* neighbour[4];
neighbour[0] = GetDipole(x,y-1); //up = 0, right = 1, down = 2, left = 3
neighbour[1] = GetDipole(x+1,y);
neighbour[2] = GetDipole(x,y+1);
neighbour[3] = GetDipole(x-1,y);
int myDirection = central->GetDirection() ;
int type = central->GetType();
if(type==0){p1 = pZero; } else{p1 = pOne; }
for(int i=0; i<4; i++){
int neighbourDirection = neighbour[i]->GetDirection();
int Ntype = neighbour[i]->GetType();
if(Ntype==1){ p2 = pOne; } else{ p2 = pZero; }
// std::cout<<"1"<<std::endl;
double Energy[4][4]=
{
{-2*k*p1*p2/r/r/r - (p1+p2)*E, -p1*E, 2*k*p1*p2/r/r/r - (p1-p2)*E, -p1*E }, // ABOVE OR BELOW DIPOLES
{-p2*E, 2*k*p1*p2/3/r/r/r, p2*E, -2*k*p1*p2/3/r/r/r },
{2*k*p1*p2/r/r/r + (p1-p2)*E, p1*E, -2*k*p1*p2/r/r/r + (p1+p2)*E, p1*E},
{-p2*E, -2*k*p1*p2/3/r/r/r, p2*E, 2*k*p1*p2/3/r/r/r}
};
// std::cout<<"2"<<std::endl;
double Energy2[4][4] =
{
{ 2*k*p1*p2/3/r/r/r -(p1+p2)*E, -p1*E, -2*k*p1*p2/3/r/r/r -(p1-p2)*E, -p1*E }, // RIGHT OR LEFT ENERGIES
{ -p2*E, -2*k*p1*p2/r/r/r, p2*E, 2*k*p1*p2/r/r/r },
{ -2*k*p1*p2/3/r/r/r +(p1-p2)*E, p1*E, 2*k*p1*p2/3/r/r/r +(p1+p2)*E, p1*E },
{-p2*E, 2*k*p1*p2/r/r/r, p2*E, -2*k*p1*p2/r/r/r }
};
//std::cout<<"3"<<std::endl;
std::cout<<"i = "<< i << "\t" << type << "\t "<<Ntype <<"\t"<< neighbourDirection<<"\t"<<myDirection<<std::endl;
std::cout<<"x = " << x << "y = "<< y <<std::endl;
if(i==0 || i==2){
energy += Energy[myDirection][neighbourDirection];
}
else if(i==1 || i==3){
energy += Energy2[myDirection][neighbourDirection];
}
}//std::cout<<"4"<<std::endl;
std::cout<<"4end clac local E"<<std::endl;
return energy;
}
Lattice::~Lattice(){
delete[] dipoleArray;
}
EDIT:: When it crashes it is always at an edge of the matrix. Usually the top row.
If my dipole is top row, it crashes when I find the 'above' neighbour.
When returning properties of this neighbour, "type" and "direction," defined in the dipole class to be 0,1 or 0,1,2,3 it returns large +ve or -ve values for the neighbour
PLEASE help as I cant see any problems :S