From d2d45e59d5437c83534ca03371a32d59af6b38e0 Mon Sep 17 00:00:00 2001
From: Pascal Engeler <engelerp@phys.ethz.ch>
Date: Thu, 16 Jan 2020 11:11:36 +0100
Subject: [PATCH] Added matrix calculation
---
include/matrix_element_calculator.hpp | 16 ++++++++++++
include/system.hpp | 36 ++++++++++++++++++++++++++-
2 files changed, 51 insertions(+), 1 deletion(-)
create mode 100644 include/matrix_element_calculator.hpp
diff --git a/include/matrix_element_calculator.hpp b/include/matrix_element_calculator.hpp
new file mode 100644
index 0000000..b8963ae
--- /dev/null
+++ b/include/matrix_element_calculator.hpp
@@ -0,0 +1,16 @@
+#ifndef MATRIX_ELEMENT_CALCULATOR_HPP_INCLUDED
+#define MATRIX_ELEMENT_CALCULATOR_HPP_INCLUDED
+
+template <value_t, params_t, vars_t, drum_t>
+class MatrixElementCalculator{
+ public:
+ MatrixElementCalculator() = default;
+ ~MatrixElementCalculator() = default;
+
+ //diagonal element at (index, index)
+ virtual value_t operator()(const size_t index, const std::vector<drum_t>& drums) const noexcept = 0;
+ //coupling element at (index1, index2)
+ virtual value_t operator()(const size_t index1, const size_t index2, const std::vector<drum_t>& drums) const noexcept = 0
+};
+
+#endif
diff --git a/include/system.hpp b/include/system.hpp
index dac90a2..2d2384f 100644
--- a/include/system.hpp
+++ b/include/system.hpp
@@ -3,7 +3,7 @@
#include <vector>
#include <iostream>
-template<typename value_t, typename drum_t, typename grabber_t, typename sysparams_t, typename force_t, typename coupler_t, typename driver_t, typename stepper_t>
+template<typename value_t, typename drum_t, typename grabber_t, typename sysparams_t, typename force_t, typename coupler_t, typename driver_t, typename stepper_t, typename matelecalc_t>
class System{
public:
System(const value_t t_end, const value_t dt, const std::vector<drum_t>& drums, const stepper_t stepper, const force_t force, const sysparams_t sysparams, const grabber_t grabber)
@@ -50,11 +50,17 @@ class System{
return grabber_.save();
}
+ std::vector<value_t> get_matrix(const matele_t& mec){
+ calculate_matrix();
+ return matrix_;
+ }
+
private:
void step_dc(const value_t dt){
sysparams_.coupler.step(dt);
sysparams_.driver.step(dt);
}
+
void push_dc(){
for(size_t i = 0; i < drums_.size()-1; ++i){
drums_[i].get_variables().V = sysparams_.driver(i);
@@ -64,6 +70,32 @@ class System{
}
}
+ //Calculates dynamical matrix for the currently pushed values
+ void calculate_matrix(const matelecalc_t& mec){
+ //clear the matrix
+ matrix_.clear();
+ //prepare memory
+ matrix_.reserve((drums_.size()-1)*(drums_.size()-1) - matrix_.capacity());
+ //initialize
+ for(size_t i = 0; i < (drums_.size()-1) * (drums_.size()-1); ++i)
+ matrix_.push_back(value_t(0.));
+
+ for(size_t i = 0; i < drums_.size()-1; ++i){
+ size_t n0_index = sysparams_.adjacency_vector[i][0];
+ size_t n1_index = sysparams_.adjacency_vector[i][1];
+ size_t n2_index = sysparams_.adjacency_vector[i][2];
+ //diagonal
+ matrix_[i][i] = mec(i, drums_);
+ //couplings
+ if(n0_index != drums_.size()-1) [[likely]] //actually a neighbour
+ matrix_[i][n0_index] = mec(i, n0_index, drums_);
+ if(n1_index != drums_.size()-1) [[likely]] //actually a neighbour
+ matrix_[i][n1_index] = mec(i, n1_index, drums_);
+ if(n2_index != drums_.size()-1) [[likely]] //actually a neighbour
+ matrix_[i][n2_index] = mec(i, n2_index, drums_);
+ }
+ }
+
std::vector<drum_t> drums_; //The last drum is NOT TO BE TOUCHED!
sysparams_t sysparams_;
grabber_t grabber_;
@@ -71,6 +103,8 @@ class System{
force_t force_;
value_t dt_, t_end_;
value_t time_;
+
+ std::vector<value_t> matrix_; //Matrix of this problem
};
#endif
--
GitLab