To receive notifications about scheduled maintenance, please subscribe to the mailing-list gitlab-operations@sympa.ethz.ch. You can subscribe to the mailing-list at https://sympa.ethz.ch

Commit 12ea4e1c authored by dkammer's avatar dkammer
Browse files

week 9

parent f9e68568
build
*~
\ No newline at end of file
# Assignment
This week we look a bit more into data structures. We work on the same code we have used over the past couple of weeks. The objective here is to find a more optimal data structure for the matrix. The direct stiffness method, much like the finite-element method, leads to a matrix that contains mostly zeros. Such matrices are referred to as Sparse Matrix. For this assignment, follow these steps:
1. read about Sparse Matrix: [https://en.wikipedia.org/wiki/Sparse_matrix](https://en.wikipedia.org/wiki/Sparse_matrix) and summarize on 1-2 slides (e.g., underlying data structure and advantages/disadvantages)
2. evaluate the changes needed in the code to use a sparse matrix, and summarize on 1-2 slides
3. make the changes and summarize problems/insights/result on slides
As a reminder: The code is provided in the folder `code`. It can be compile via `cmake`. Follow these steps:
1. `cd code`
2. `mkdir build`
3. `cd build`
4. `ccmake ..`
5. configure and generate
6. `make`
cmake_minimum_required(VERSION 3.1)
project(tehpc)
enable_language(CXX)
# flags for compilation (options)
set(CMAKE_CXX_FLAGS "-Wall -Wextra -pedantic -std=c++11 " CACHE STRING "" FORCE)
set(CMAKE_CXX_FLAGS_RELEASE "-O3 -DNDEBUG " CACHE STRING "" FORCE )
# alternative to add_library (but with option)
set(BUILD_SHARED_LIBS ON CACHE BOOL "Build shared libraries.")
mark_as_advanced(BUILD_SHARED_LIBS)
# build type
set(CMAKE_BUILD_TYPE "${CMAKE_BUILD_TYPE}" CACHE STRING
"Build options: None Debug Release RelWithDebInfo MinSizeRel."
FORCE )
# include directories
set(TEHPC_INCLUDE_DIRS
"${CMAKE_CURRENT_BINARY_DIR}/src"
"${CMAKE_CURRENT_SOURCE_DIR}/src"
)
include_directories(${TEHPC_INCLUDE_DIRS})
find_package(Armadillo REQUIRED)
include_directories(${ARMADILLO_INCLUDE_DIRS})
# source directory
add_subdirectory(src)
# ------------------------------------------------------
# examples
# define add_simulation function
function(add_simulation SIM_EXE)
add_executable(${SIM_EXE} ${SIM_EXE}.cpp)
target_link_libraries(${SIM_EXE} tehpc)
endfunction()
option(TEHPC_EXAMPLES "examples" OFF)
if(TEHPC_EXAMPLES)
add_subdirectory(examples)
endif()
# ------------------------------------------------------
# test
option(TEHPC_TEST "tests" OFF)
if(${TEHPC_TEST})
enable_testing()
add_subdirectory(tests)
endif()
This diff is collapsed.
add_simulation(basic_example)
configure_file(basic_example_1.inp basic_example_1.inp COPYONLY)
configure_file(basic_example_1.run basic_example_1.run COPYONLY)
#include <iostream> // print to console (cout)
//#include <cmath> // math operations
//#include <fstream> // file streaming for input file
#include "Model.hpp"
int main(int argc, char *argv[])
{
std::string fname;
std::string odir;
if(argc<3) {
std::cerr << "Not enough arguments:"
<< " ./basic_example <sim-name> <output-directory>" << std::endl;
return 0;
}
else {
fname = argv[1];
odir = argv[2];
}
Model my_model;
std::cout << "read input" << std::endl;
my_model.read_input(fname+".inp");
my_model.assemble(); // Get unconstrained stiffness matrix + load vector
my_model.apply_BC();
std::cout << "solve" << std::endl;
my_model.solve(); // Solve linear system and get displacement vector
my_model.compute_reaction(); // Calculate reaction
my_model.print_results(odir,fname); // print results
return 0;
}
2 # Nel (Number of elements)
2 # Nen (Number of nodes per element)
3 # Nnp (Number of nodes)
2 # n (Number of dofs / degrees of freedom per node)
1 0; 0 0; 1 1 # nodes (node coordinates (x,y))
1 3; 2 3 # elements (Nodes connecting each element)
1 1 # E (Young Modulus of all elements)
1 1 # A (Cross-section area of all elements)
0 0; 0 0; 1 1 # pli (Prescribed Load Index)
0 0; 0 0; 10 0 # pl (Prescribed Load)
1 1; 1 1; 0 0 # pdi (Prescribed Displacement Index)
0 0; 0 0; 0 0 # pd (Prescribed Displacement)
odir=~/dev/output
name=basic_example_1
./basic_example $name $odir 2>&1 | tee $name.progress
cp $name.inp $odir/$name.inp
cp $name.run $odir/$name.run
cp $name.progress $odir/$name.progress
set(TEHPC_HEADER
Model.hpp
#Mesh.hpp
)
set(TEHPC_SRC
Model.cpp
#Mesh.cpp
)
set(TEHPC_DEPEND_FILES ${TEHPC_SRC} ${TEHPC_HEADER} PARENT_SCOPE)
add_library(tehpc ${TEHPC_SRC})
# if linking to other libraries, e.g. MPI
target_link_libraries(tehpc ${ARMADILLO_LIBRARIES})
#include <iostream> // print to console (cout)
#include <cmath> // math operations
#include "Model.hpp"
#define LEN 256 // buffer length for input reading
// templated function to read any armadillo object from line string.
// (armadillo constructor from string -> ok for small matrices but not efficient).
template <typename T>
void Model::read_line(T * el, std::ifstream &fin){
char line[LEN]; // Buffer for storing a line from the file
fin.getline(line,LEN,'#'); // Load line from file until '#' character is reached
*el = T(line); // Convert 'line' string to T type (mat, imat, vec, etc)
fin.getline(line,LEN); // Move the file stream pointer to the next line.
}
// read integer from line.
void Model::read_line(int * el, std::ifstream &fin){
char line[LEN];
fin.getline(line,LEN,'#');
*el = atoi(line);
fin.getline(line,LEN);
}
// INPUT FUNCTION : Reads input from file. Takes file name and returns input struct.
// ---------------------------------------------------------------------------------
void Model::read_input(const std::string & file){
std::ifstream fin(file); // Loads file stream
read_line(&Nel, fin);
read_line(&Nen, fin);
read_line(&Nnp, fin);
read_line(&n, fin);
read_line<arma::mat>(&nodes, fin);
read_line<arma::imat>(&elements, fin);
read_line<arma::vec>(&E, fin);
read_line<arma::vec>(&A, fin);
read_line<arma::imat>(&pli, fin);
read_line<arma::mat>(&pl, fin);
read_line<arma::imat>(&pdi, fin);
read_line<arma::mat>(&pd, fin);
}
void Model::local_stiffness_element(arma::mat & k_e,
int el) {
double dx = nodes(elements(el,1)-1,0) - nodes(elements(el,0)-1,0);
double dy = nodes(elements(el,1)-1,1) - nodes(elements(el,0)-1,1);
double L = sqrt(dx*dx + dy*dy);
double cx = dx / L; // Sin of element angle with respect to global coordinates
double cy = dy / L; // Cosin of element angle with respect to global coordinates
double k = (A(el)*E(el))/L;
// local K
arma::mat k_e_local( {{ k, 0, -k, 0},
{ 0, 0, 0, 0},
{-k, 0, k, 0},
{ 0, 0, 0, 0}});
// Rotation matrix from local angle to global
arma::mat R ({{ cx, cy, 0, 0},
{-cy, cx, 0, 0},
{ 0, 0, cx, cy},
{ 0, 0, -cy, cx}});
arma::mat rhs = k_e_local * R; // Store right-hand-side
k_e = R.t()*rhs; // Rotated local K of element i (still local indeces)
// R.t() -> Transpose of R. For orthogonal matrices: R^-1 = R^T
}
//---------------------------------------------------------------------------------------
// ASSEMBLY FUNCTION
// Assemble stiffness matrix from input struct containing all parsed input values from input file.
// -----------------------------------------------------------------------------------
void Model::assemble() {
arma::imat LM(Nen*n,Nel); // Gather arma::matrix --> Converts local indeces of dofs to global
for (int i=0; i<Nel; ++i) {
for (int j=0; j<Nen; ++j) {
for (int k=0; k<n; ++k) {
LM(2*j+k,i) = k + n * elements(i,j) - 2;
}
}
}
// Floating stiffness matrix (No boundary conditions / )
arma::mat* p_k_f = & K_f;
arma::mat K_f(Nnp*n,Nnp*n,arma::fill::zeros);
// Assemble floating K
for (int i=0;i<Nel;++i) {
arma::mat k_e;
this->local_stiffness_element(k_e, i);
arma::mat K_e(Nnp*n,Nnp*n,arma::fill::zeros);
for (int j=0;j<Nen*n;++j) {
for (int k=0;k<Nen*n;++k) {
K_e(LM(j,i),LM(k,i)) = k_e(j,k); // Re-arrange local k_e indices into global K
}
}
K_f = K_f + K_e; // Add the rotated and globally adjusted dofs of element i to global K
}
/* Swap exchanges the content the two objects with each other.
* We could also write `model->K = K` which would copy the matrix named `K`
* into the member `K` of the `model` object, however this would copy the matrix.
* Since we do not need the matrix any longer, we can do this.
*/
(*p_k_f).swap(K_f);
}
// applies boundary conditions on displacements
void Model::apply_BC (){
arma::mat* p_k = & K;
arma::mat* p_f = & f;
arma::mat* p_pd = & pd;
arma::mat K (K_f); // Start from floating structure
arma::vec f (pl.as_row().as_col()); // Initialize load
arma::uvec pdl = find(pdi.as_row()); // contains 1 for prescribed, 0 for free (displacement dofs)
arma::vec pd ((*p_pd).as_row().as_col()); // prescribed displacements.
for (int i=0;(unsigned)i< pdl.n_elem;++i){
f(pdl(i)) = pd(pdl(i));
K.row(pdl(i)) = arma::rowvec (Nnp*n,arma::fill::zeros);
K(pdl(i),pdl(i)) = 1;
}
(*p_k).swap(K);
(*p_f).swap(f);
}
// REACTION FUNCTION : Computes reaction to load. Takes floating stiffness matrix
// + displacements + input parameters (prescriptions) and returns the reaction vector.
void Model::compute_reaction (){
arma::vec* p_r = & r;
arma::vec r = K_f * d;
arma::uvec pdl = find(pli.as_row());
for (int i=0;(unsigned)i<pdl.n_elem;++i){
r(pdl(i)) = 0;
}
(*p_r).swap(r);
}
void Model::solve()
{
*(& d) = arma::solve (K,f); // Solve linear system and get displacement vector
}
void Model::print_results(const std::string & path,
const std::string & fname)
{
std::string path_to_file = path + "/" + fname + ".out";
std::ofstream ofile(path_to_file, std::ios::out);
ofile << "K_f:\n" << K_f << "\n";
ofile << "K:\n" << K << "\n";
ofile << "d:\n" << d << "\n";
ofile << "r:\n" << r << "\n";
ofile.close();
}
#ifndef MODELHEADERDEF
#define MODELHEADERDEF
#include <fstream> // file streaming for input file
#include <armadillo> // armadillo library
class Model
{
public:
template <typename T>
void read_line(T*, std::ifstream&);
void read_line(int*, std::ifstream&);
void read_input(const std::string&);
void assemble ();
void apply_BC ();
void solve ();
void compute_reaction ();
void print_results(const std::string & path,
const std::string & fname);
protected:
void local_stiffness_element(arma::mat & k_e, int el);
private:
int Nel,Nen,Nnp,n; // Number of elements, nodes-per-element, nodes, dofs per node
arma::mat nodes; // Coordinates of all nodes
arma::imat elements; // Node index pair from every element
arma::vec E, A; // Young modulus and cross section of every element
arma::imat pli, pdi; // constrained dofs matrix (1 -> constrained / 0 -> free)
arma::mat pl, pd; // load and displacement boundary conditions for every dof
// Struct to store all vectors and matrices of the 2d truss stiffness model.
arma::mat K,K_f; // Final stiffness matrix + unconstrained stiffness matrix (floating).
arma::vec f,d,r; // Load + displacement + reaction.
};
#endif
# define a function to implement test
function(tehpc_add_test TEST_NAME)
add_simulation(${TEST_NAME})
add_test(${TEST_NAME} ${TEST_NAME})
endfunction()
# dummy test
#tehpc_add_test(dummy.test)
#include <iostream>
int main()
{
std::cout << "start dummy test" << std::endl;
/*
if (some-condition) {
std::cerr << "test failed" << std::endl;
return 1; // failure
}
*/
// all tests passed: success
std::cout << "end dummy test" << std::endl;
return 0; // success
}
build
*~
\ No newline at end of file
# Assignment
This week we look a bit more into data structures. We work on the same code we have used over the past couple of weeks. The objective here is to learn about iterators and use them to loop over built-in data structures. For this assignment, follow these steps:
1. read about iterators:
a. [https://www.geeksforgeeks.org/introduction-iterators-c/](https://www.geeksforgeeks.org/introduction-iterators-c/)
b. [https://www.geeksforgeeks.org/difference-between-iterators-and-pointers-in-c-c-with-examples/](https://www.geeksforgeeks.org/difference-between-iterators-and-pointers-in-c-c-with-examples/)
c. [https://www.cplusplus.com/reference/iterator/](https://www.cplusplus.com/reference/iterator/)
2. summarize what you learnt on 1-2 slides
3. evaluate how you could use iterators in the `velverlet` function and summarize on 1-2 slides
4. make the changes and summarize problems/insights/results on slides (keep a copie of for the moment)
5. measure the time for the `velverlet` and compare to previous code. Helpful code to measure time:
a. `#include <chrono> // to time function execution`
b. `auto t1 = std::chrono::high_resolution_clock::now(); // get current time`
c. `auto dt = std::chrono::duration_cast<std::chrono::microseconds>(t2 - t1).count() // time passed between t1 and t2`
As a reminder: The code is provided in the folder `code`. It can be compile via `cmake`. Follow these steps:
1. `cd code`
2. `mkdir build`
3. `cd build`
4. `ccmake ..`
5. configure and generate
6. `make`
build
\ No newline at end of file
cmake_minimum_required(VERSION 3.1)
project(tehpc)
enable_language(CXX)
# flags for compilation (options)
set(CMAKE_CXX_FLAGS "-Wall -Wextra -pedantic -std=c++11 " CACHE STRING "" FORCE)
set(CMAKE_CXX_FLAGS_RELEASE "-O3 -DNDEBUG " CACHE STRING "" FORCE )
# alternative to add_library (but with option)
set(BUILD_SHARED_LIBS ON CACHE BOOL "Build shared libraries.")
mark_as_advanced(BUILD_SHARED_LIBS)
# build type
set(CMAKE_BUILD_TYPE "${CMAKE_BUILD_TYPE}" CACHE STRING
"Build options: None Debug Release RelWithDebInfo MinSizeRel."
FORCE )
# include directories
set(TEHPC_INCLUDE_DIRS
"${CMAKE_CURRENT_BINARY_DIR}/src"
"${CMAKE_CURRENT_SOURCE_DIR}/src"
)
include_directories(${TEHPC_INCLUDE_DIRS})
find_package(Armadillo REQUIRED)
include_directories(${ARMADILLO_INCLUDE_DIRS})
# source directory
add_subdirectory(src)
# ------------------------------------------------------
# examples
# define add_simulation function
function(add_simulation SIM_EXE)
add_executable(${SIM_EXE} ${SIM_EXE}.cpp)
target_link_libraries(${SIM_EXE} tehpc)
endfunction()
option(TEHPC_EXAMPLES "examples" OFF)
if(TEHPC_EXAMPLES)
add_subdirectory(examples)
endif()
# ------------------------------------------------------
# test
option(TEHPC_TEST "tests" OFF)
if(${TEHPC_TEST})
enable_testing()
add_subdirectory(tests)
endif()
This diff is collapsed.
add_simulation(basic_example)
configure_file(basic_example_1.inp basic_example_1.inp COPYONLY)
configure_file(basic_example_1.rest basic_example_1.rest COPYONLY)
configure_file(basic_example_1.run basic_example_1.run COPYONLY)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment