Generalized Potential, general refactoring, modified tests

src/Basis/Base.cpp

@@ -35,3 +35,24 @@ std::vector<ContinuousBase> Base::getContinuous() {
 std::vector<DiscreteBase> Base::getDiscrete() {
 	return this->discrete;
 }
+
+std::ostream& operator<<(std::ostream& stream, Base& base) {
+
+	// Print continuous dimension values (if present)
+	if (base.getContinuous().size() > 0)
+		for (int i = 0; i < base.getContinuous().size(); i++) {
+			for (int coord_counter = 0; coord_counter < base.getContinuous().at(i).getCoords().size(); coord_counter++)
+				stream << base.getContinuous().at(i).getCoords().at(coord_counter) << "; ";
+		}
+
+	stream << "\n\n";
+
+	// Print discrete dimension values (if present)
+	if (base.getDiscrete().size() > 0)
+		for (int i = 0; i < base.getContinuous().size(); i++) {
+			for (int coord_counter = 0; coord_counter < base.getDiscrete().at(i).getCoords().size(); coord_counter++)
+				stream << base.getDiscrete().at(i).getCoords().at(coord_counter) << "; ";
+		}
+
+    return stream;
+ }

src/Basis/Base.h

@@ -19,6 +19,7 @@ class Base
 	int getDim();
 	std::vector<ContinuousBase> getContinuous();
 	std::vector<DiscreteBase> getDiscrete();
+	friend std::ostream& operator<< (std::ostream& stream, Base& base);
 
 private:
 	std::vector< DiscreteBase > discrete;

src/Basis/BasisManager.cpp

@@ -21,8 +21,8 @@ void BasisManager::addBase(Base b) {
 	this-> bases.push_back(b);
 
 	// If it's there's just this one in the vector, then it's automatically selected
-		if (this->bases.size() == 1)
-		this->selectBase(b);
+	if (this->bases.size() == 1)
+	this->selectBase(b);
 }
 
 std::vector<Base> BasisManager::getBasisList() {
@@ -57,20 +57,19 @@ Base BasisManager::Builder::build(Base::basePreset b, int dimension, double mesh
     //TODO: Eventually add controls...
 
     switch(b){
-
-        case Base::basePreset::Custom: throw std::invalid_argument("Custom basis not meaningful with parameters!");
-            break;
         case Base::basePreset::Cartesian:
             std::cout<< "Building Cartesian basis in " << dimension << "dimensions, with nbox = " << nbox << ", mesh = " << mesh << std::endl;
-            for(int i = 0; i < dimension; i++)
-            {
+            for(int i = 0; i < dimension; i++) {
                 this->addContinuous(mesh, nbox);
             }
             break;
+
         // todo: find criterium criterium to define Lmax given only r. Harmonic oscillator r/l relation?
-        case Base::basePreset::Cylindrical: std::invalid_argument("Wrong parameters for Cylindrical basis");
+        case Base::basePreset::Cylindrical: throw std::invalid_argument("Wrong parameters for Cylindrical basis");
+            break;
+        case Base::basePreset::Spherical: throw std::invalid_argument("Wrong parameters for Spherical basis");
             break;
-        case Base::basePreset::Spherical: std::invalid_argument("Wrong parameters for Spherical basis");
+        case Base::basePreset::Custom: throw std::invalid_argument("Custom basis not meaningful with parameters!");
             break;
         default: throw std::invalid_argument("Wrong basis type or initialization meaningless!");
             break;
@@ -92,16 +91,14 @@ Base BasisManager::Builder::build(SphericalInitializer ini) {
 }
 
 Base BasisManager::Builder::build(ContinuousInitializer ini) {
+
 	if(ini.mesh <= 0 ) std::invalid_argument("mesh < 0 does not have sense");
 	if(ini.end <= ini.mesh ) std::invalid_argument("xmax < mesh does not have sense");
-
-	std::cout<< "Building Cartesian basis between " << ini.start << ", " << ini.end << "  mesh = " << ini.mesh << std::endl;
-
 	this->addContinuous(ini.start,ini.end,ini.mesh);
 
+	std::cout<< "Building Cartesian basis between " << ini.start << ", " << ini.end << "  mesh = " << ini.mesh << std::endl;
 	return Base(Base::basePreset::Spherical, 3, c_base, d_base);
 }
-
 // --- End Factory --- //
 
 

src/Basis/BasisManager.h

@@ -1,3 +1,6 @@
+#ifndef BASISMANAGER_H
+#define BASISMANAGER_H
+
 #include <Base.h>
 #include <Initializer.h>
 
@@ -38,3 +41,4 @@ class BasisManager
 	static BasisManager* instance;
 	BasisManager() {}
 };
+#endif

src/Potential/Potential.cpp

@@ -1,54 +1,101 @@
 #include "Potential.h"
 
-Potential::Potential(Base base, std::string type, double k, double width, double height)
+Potential::Potential(Base base, PotentialType type, double k, double width, double height, bool separable)
 {
-    this->base     = base; 
-    this->k        = k;
-    this->width    = width;
-    this->height   = height;
-    this->type     = type;
+    this->base      = base; 
+    this->k         = k;
+    this->width     = width;
+    this->height    = height;
+    this->type      = type;
+    this->separable = separable;
+    this->v         = std::vector<double>();
 
-    // The following is hardcoded to work with a base having one continuous dimension. 
-    this->v        = std::vector<double>(base.getContinuous().at(0).getCoords().size());
+    if (!this->separable) {
+        switch(type) {
+            case BOX_POTENTIAL:
+                this->box_potential();
+                break;
+            case HARMONIC_OSCILLATOR:
+                this->ho_potential();
+                break;
+            case FINITE_WELL_POTENTIAL:
+                this->finite_well_potential();
+                break;
+            default:
+                throw std::invalid_argument("Wrong potential type or initialization meaningless!");
+        }
+    }
+    // If the potential is separable, then n potentials (one for each base's representation)
+    else {
+        this->separated_potentials = std::vector<Potential>();
 
-    if(type.compare("box potential") == 0 || type.compare("box") == 0 || type.compare("0") == 0)
-        this->box_potential();
+        for (int i = 0; i < base.getContinuous().size(); i++) {
+            // Extract continuous representation from the main base
+            std::vector<ContinuousBase> c_base = std::vector<ContinuousBase>(1);
+            c_base.push_back(base.getContinuous().at(i));
 
-    else if(type.compare("harmonic oscillator") == 0 || type.compare("ho") == 0 ||  type.compare("1") == 0)
-        this->ho_potential();
+            // Create new base object with the extracted representation
+            Base monodimensional = Base(Base::basePreset::Custom, 1, c_base, std::vector<DiscreteBase>());
 
-    else if(type.compare("finite well potential") == 0 || type.compare("well") == 0 || type.compare("2") == 0 )
-        this->finite_well_potential();
+            // Create new potential with the new base
+            Potential separated_potential = Potential(monodimensional, this->type, this->k, this->width, this->height, false);
 
-    else {
-        throw std::invalid_argument("Wrong potential type ("+type+") or initialization meaningless!");
+            // Add the new Potential to the separated potenial vector associated to the main potential
+            this->separated_potentials.push_back(separated_potential);
+        }
+
+        for (int i = 0; i < base.getDiscrete().size(); i++) {
+            // Extract discrete representation from the main base
+            std::vector<DiscreteBase> d_base = std::vector<DiscreteBase>(1);
+            d_base.push_back(base.getDiscrete().at(i));
+
+            // Create new base object with the extracted representation
+            Base monodimensional = Base(Base::basePreset::Custom, 1, std::vector<ContinuousBase>(), d_base);
+
+            // Create new potential with the new base
+            Potential separated_potential = Potential(monodimensional, this->type, this->k, this->width, this->height, false);
+
+            // Add the new Potential to the separated potenial vector associated to the main potential
+            this->separated_potentials.push_back(separated_potential);
+        }
     }
 }
 
 void Potential::ho_potential()
 {
- /*
- for(std::vector<int>::size_type i = 0; i < x.size(); i++)
-        this->v[i] = this->x[i] * this->x[i] * this->k;
-    */
+    std::vector<double> x = this->getCoordsFromBase();
+
+    for(std::vector<int>::size_type i = 0; i < x.size(); i++)
+            this->v[i] = x[i] * x[i] * this->k;
 }
 
 void Potential::box_potential()
 {
- //   std::fill(this->v.begin(), this->v.end(), 0.0);
+    std::fill(this->v.begin(), this->v.end(), 0.0);
 }
 
 void Potential::finite_well_potential()
 {
-/*
+    std::vector<double> x = this->getCoordsFromBase();
+
     for(std::vector<int>::size_type i = 0; i < x.size(); i++)
-        this->v[i] = (this->x[i] > -this->width/2.0 && this->x[i] < this->width/2.0) ? 0.0 : this->height;
-*/
+        this->v[i] = (x[i] > -this->width/2.0 && x[i] < this->width/2.0) ? 0.0 : this->height;
 }
 
-
-
 std::vector<double> Potential::getValues()
 {
     return this->v;
 }
+
+std::vector<double> Potential::getCoordsFromBase() 
+{
+    if (this->base.getContinuous().size() == 1)
+        return this->base.getContinuous().at(0).getCoords();
+
+    else if (this->base.getDiscrete().size() == 1) {
+        // tricky conversion taking each std::vector<int> value and returning a final std::vector<double> 
+        std::vector<int> original_coords = this->base.getDiscrete().at(0).getCoords();
+        return std::vector<double>(original_coords.begin(), original_coords.end());
+    }
+}
+

src/Potential/Potential.h

@@ -6,6 +6,7 @@
 #include <string>
 #include <algorithm>
 #include <stdexcept>
+#include <stdbool.h>
 #include "../Basis/Base.h"
 
 /*! Class Potential contains the potential used in the Schroedinger equation.
@@ -28,41 +29,55 @@
  */
 
 class Potential {
-private:
-    Base base;
-    std::vector<double> v;
-    std::string type;
-
-    double k;
-    double width;
-    double height;
-
-    void ho_potential();
-    void box_potential();
-    void finite_well_potential();
 
 public:
-    Potential(Base, std::string, double, double, double);
+    enum PotentialType {
+        BOX_POTENTIAL,
+        HARMONIC_OSCILLATOR,
+        FINITE_WELL_POTENTIAL,   
+    };
+
+    Potential(Base, PotentialType, double, double, double, bool);
     std::vector<double> getValues();
+    std::vector<double> getCoordsFromBase();
     Base getBase();
 
     class Builder{
         private:
             Base base;
-            std::string type     = "box";
+            PotentialType type   = PotentialType::BOX_POTENTIAL;
             double k             = 0.5;
             double width         = 5.0;
             double height        = 10.0;
+            bool separable       = false;
+            std::vector<Potential> separated_potentials;
 
         public:
             Builder(Base b);
             Builder setK(double k_new);
             Builder setWidth(double width_new);
             Builder setHeight(double height_new);
-            Builder setType(std::string type);
+            Builder setType(PotentialType type);
             Builder setBase(Base b);
+            Builder setSeparable(bool separable);
             Potential build();
     };
+
+private:
+    Base base;
+    std::vector<double> v;
+    PotentialType type;
+
+    double k;
+    double width;
+    double height;
+    bool separable;
+    std::vector<Potential> separated_potentials;
+
+    void ho_potential();
+    void box_potential();
+    void finite_well_potential();
+
 };
 
 #endif

src/Potential/PotentialBuilder.cpp

@@ -26,18 +26,25 @@ Potential::Builder Potential::Builder::setHeight(double height_new)
     return *this;
 }
 
-Potential::Builder Potential::Builder::setType(std::string type)
-{
-    if (!type.empty()) {
+Potential::Builder Potential::Builder::setType(PotentialType type)
+{   
+    // When you add new potential types in the enum, change this check!
+    if (type >= PotentialType::BOX_POTENTIAL && type <= PotentialType::HARMONIC_OSCILLATOR) {
         this->type = type;
         return *this;
     }
-    else throw std::invalid_argument("Empty type given as parameter.");
+    else throw std::invalid_argument("Wrong type given as parameter.");
+}
+
+Potential::Builder Potential::Builder::setSeparable(bool separable)
+{
+    this->separable = separable;
+
 }
 
 Potential Potential::Builder::build(){
     try {
-        return Potential(this->base,this->type,this->k,this->width,this->height);
+        return Potential(this->base,this->type,this->k,this->width,this->height, this->separable);
     }
     catch(const std::invalid_argument& e){
         throw;

src/main.cpp

@@ -30,7 +30,7 @@ int main(int argc, char **argv) {
     builder.addContinuous(start_1, end_2, mesh_2);
     builder.addContinuous(start_2, end_2, mesh_2);
 
-	// Getting a base object, building the previous specified in the basis builder.
+	// Getting a base object, building with the previous specified parameters in the basis builder.
 	// (Specify that the base is Cartesian, passing also the dimension)
 	base = builder.build(Base::basePreset::Cartesian, dimension);
 
@@ -41,21 +41,7 @@ int main(int argc, char **argv) {
 	basis = manager->getBasisList();
 
     // Print basis values for each dimension
-	for (int dimension_counter = 0; dimension_counter < base.getDim(); dimension_counter++) {
-        std::cout << "Dimension: " << dimension_counter+1 << std::endl;
-
-        // Print continuous dimension values (if present)
-        if (base.getContinuous().size() > 0)
-            for (int coord_counter = 0; coord_counter < base.getContinuous().at(dimension_counter).getCoords().size(); coord_counter++)
-                std::cout << base.getContinuous().at(dimension_counter).getCoords().at(coord_counter) << "; ";
-
-        // Print discrete dimension values (if present)
-        if (base.getDiscrete().size() > 0)
-            for (int coord_counter = 0; coord_counter < base.getDiscrete().at(dimension_counter).getCoords().size(); coord_counter++)
-                std::cout << base.getDiscrete().at(dimension_counter).getCoords().at(coord_counter) << "; ";
-
-        std::cout << std::endl << std::endl;
-    }
+	std::cout << base;
 
     return 0;
 }