計算屋さんの日常

2014/01/29

高次の形状関数

有限要素法の高次形状関数の実装は, 出来そうだ.
c++ならばtemplateをつかって再帰的に計算できる.

参考にしました.
http://www.colorado.edu/engineering/cas/courses.d/IFEM.d/IFEM.Ch18.d/IFEM.Ch18.pdf

高次のガウス積分の実装を今度は考えなくては.

2014/01/27

static polymorphismとgcc

はてgccのバージョン4.8.2にしたら静的多様性を使ったコードがコンパイル出来なくなったぞい. なぜじゃ?　4.7.3でもでるな　???? なんかミスったか ?

c/c++, fortranとかは普通に作れるんだけど, meta化すると急にできなくなる.

2014/01/25

とある先進的な解析ソフトを見つけてコンパイルしようとしたら, どうもconfigすら通らない. いろいろ探ってみたのだけれども理由がわからず, mailing listに質問してみたら

"ごめん, ホームページ更新してなくて. 依存情報かなり古いんだと言われた"

おっさんは, その程度でめげる人ではないので,

gmp-5.1.3
mpfr-3.1.2
mpc-1.0.2
gcc-4.8.2
boost-1.55.0
openmpi-1.6.5

をつくった.

あとpetsc-3.4.3を作れば問題ないはず.

2014/01/24

FLTK

FLTK -fPICオプション付けないとコンパイル通らないな. FLTK-1.3ね

2014/01/22

c++ テンプレートメタプログラミング　3 静的多様性を用いた行列の足し算

いくつでも行列足し算できる奴も出来た. gcc-4.8.0で動作確認, 4.3.4ではコンパイル出来なかった.

元ネタ
http://itee.uq.edu.au/~conrad/misc/sanderson_templates_lecture_uqcomp7305.pdf

#include<iostream>

template< typename T1, typename T2 >
class Glue;

template< typename derived>
struct Base {
   const derived & get_ref() const{
       return static_cast<const derived&>(*this);
   }
};

class Matrix : public Base < Matrix >{
public:
   int rows;
   int cols;
   double *data;

   Matrix(){}
   Matrix( int in_rows, int in_cols ){
       set_size( in_rows, in_cols );
   }
   void set_size( int in_rows, int in_cols ){
       rows = in_rows;
       cols = in_cols;
       data = new double [rows*cols];
       for( int i=0; i<rows*cols; i++ ) data[i] = 1.0;
   }

   Matrix( const Matrix & X );
   const Matrix& operator=( const Matrix & X );

   template< typename T1, typename T2 >
   Matrix( const Glue<T1, T2> & X );

   template< typename T1, typename T2 >
   const Matrix & operator=( const Glue<T1, T2> & X );
};

template< typename T1, typename T2 >
class Glue : public Base< Glue<T1,T2> >{
public:
   const T1 & A;
   const T2 & B;

   Glue( const T1 & in_A, const T2 & in_B ): A(in_A), B(in_B) {

   }
};

template< typename T1, typename T2 >
inline const Glue< T1, T2 >
operator+( const Base<T1>&A, const Base<T2> & B ){
   return Glue< T1, T2 >( A.get_ref(), B.get_ref() );
}

template< typename T1 >
struct depth_lhs{
   static const int num = 0;
};

template< typename T1, typename T2 >
struct depth_lhs< Glue<T1, T2 > >{
   static const int num = 1 + depth_lhs<T1>::num;
};

template< typename T1 >
struct mat_ptrs{
   static const int num = 0;

   inline static void
   get_ptrs( const Matrix** ptrs, const T1 & X ){
       {
           ptrs[0] = reinterpret_cast<const Matrix*>(&X);
       }
   }
};

template< typename T1, typename T2 >
struct mat_ptrs< Glue<T1, T2> >
{
   static const int num = 1 + mat_ptrs<T1>::num;

   inline static void
   get_ptrs( const Matrix** in_ptrs, const Glue<T1,T2>& X ){
       mat_ptrs<T1>::get_ptrs(in_ptrs, X.A );

       in_ptrs[num] = reinterpret_cast< const Matrix*>( & X.B );
   }
};

template< typename T1, typename T2 >
const Matrix & Matrix::operator=( const Glue<T1,T2>& X ){
   int i,j;
   double sum;
   int N = 1 + depth_lhs < Glue<T1,T2> >::num;
   const Matrix* ptrs[N];
   mat_ptrs< Glue<T1,T2> >::get_ptrs(ptrs,X);

   int r = ptrs[0]->rows;
   int c = ptrs[0]->cols;

   set_size( r, c );

   for( j=0; j<r*c; ++j ){
       double sum = ptrs[0]->data[j];
       for( i=1; i<N; ++i ){
           sum += ptrs[i]->data[j];
       }
       data[j] = sum;
   }
   return *this;
}

int main(){
   Matrix A(2,2);
   Matrix B(2,2);
   Matrix C(2,2);
   Matrix D(2,2);
   Matrix X(2,2);

   X = A+B+C+D;

   for( int i=0; i<4; i++ ) std::cout << X.data[i] << std::endl;
return 0;}

c++ テンプレートメタプログラミング 2 行列足し算

3つの行列を足し算する版できた.

元ネタ
http://itee.uq.edu.au/~conrad/misc/sanderson_templates_lecture_uqcomp7305.pdf

#include<iostream>

template < typename T1, typename T2 >
class Glue;

class Matrix {
    public:
    int rows;
    int cols;
    double *data;

    Matrix(){}
    Matrix( int in_row, int in_cols ){
        set_size( in_row, in_cols );
    }
    void set_size( int in_rows, int in_cols ){
        data = new double [ in_rows*in_cols ];
        rows = in_rows;
        cols = in_cols;
        for( int i=0; i<rows*cols; i++ ) data[i]=1.0;
    }

    Matrix( const Matrix & X );
    const Matrix operator=( const Matrix& X );

    Matrix( const Glue<Matrix,Matrix> & X );
    const Matrix & operator=( const Glue<Matrix,Matrix> & X );

    Matrix( const Glue< Glue<Matrix,Matrix>, Matrix > & X );
    const Matrix & operator=( const Glue< Glue<Matrix,Matrix>, Matrix > & X );

};

template < typename T1, typename T2 >
class Glue {
public:
    const T1 & A;
    const T2 & B;
    Glue( const T1 & in_A, const T2 & in_B ): A(in_A), B(in_B){}
};

Matrix::Matrix( const Glue<Matrix,Matrix> & X ){
    operator=(X);
}

const Matrix & Matrix::operator=( const Glue< Glue<Matrix,Matrix>, Matrix> & X ){
    const Matrix & A = X.A.A;
    const Matrix & B = X.A.B;
    const Matrix & C = X.B;

    set_size( A.rows, A.cols );

    for( int i=0; i<A.rows * A.cols; i++ ){
        data[i] = A.data[i] + B.data[i] + C.data[i];
    }
    return *this;
}

const Matrix & Matrix::operator=( const Glue<Matrix,Matrix> & X ){
    const Matrix & A = X.A;
    const Matrix & B = X.B;

    set_size( A.rows, A.cols );

    for( int i=0; i<A.rows * A.cols; i++ ){
        data[i] = A.data[i] + B.data[i];
    }
    return *this;
}

const Glue<Matrix,Matrix> operator+( const Matrix & A, const Matrix & B ){
    return Glue<Matrix, Matrix>(A,B);
}

const Glue< Glue<Matrix,Matrix> ,Matrix> operator+( const Glue<Matrix, Matrix> & P, const Matrix & Q ){
    return Glue< Glue<Matrix,Matrix>, Matrix >(P,Q);
}

int main(){
    Matrix A1(3,4);
    Matrix B1(3,4);
    Matrix C1(3,4);
    Matrix X(3,4);

    X = A1 + B1;
    X = A1 + B1 + C1;

    for( int i=0; i<X.rows*X.cols; i++ ) std::cout << X.data[i] << std::endl;

return 0;}

c++でテンプレートメタプログラミング

自作の有限要素ソルバーを, テンプレートメタプログラミングを使用して改造しようとしているがどこがどう早くなるのかよくわかっていない.

アセンブリにして吐いてみるとたしかに短くなってわいるのだけど, メタ化した部分か消えてるのかよくわからない.

アセンブリ言語読めるようになるしかないのか.

コンパイル時に計算するらしいけど, 行列式の解などは計算されると困る. 並列計算機を使ってもものすごく時間のかかる箇所なので, コンパイラーなどにやらせたら一生終わらんだろう.

そういや昔こんなの書いたね
http://a-daily-life-in-the-office.blogspot.jp/2013/02/c.html

2014/01/21

低血圧

109/69 mmHg 心拍数 66だった.

低血圧だよ.

2014/01/20

バグ

やっとバグを見つけた. 再計算投げたし, もう寝よう

2014/01/19

研究室

今所属している研究室は, 3月で教員が定年退官の為閉鎖になる. 同じ部屋に他分野の研究津市があり, 2つの分野の学生が1部屋にいる.

私ともう一人の同僚が, この研究室からの最後の博士後期課程の卒業生. そして, 後三年後に部屋をシェアしている研究室から最初の博士が出る予定.

2014/01/17

c++の二重クラス

c++のレベルがまた上がった, というかこの程度のことが知らなかったのが恥ずかしい.

二重クラスの定義

class MainClass{
public:
    MainClass();
    ~MainClass();
    class SubClass {
        public:
            SubClass();
            ~SubClass(){}
    };
};

MainClass::MainClass(){
    SubClass subclass;
}

MainClass::~MainClass(){
}

MainClass::SubClass::SubClass(){

}

int main(){
    MainClass mainclass;
    MainClass::SubClass subclass;
return 0;}

なるほど大規模なソフトでよく見る奴はこうしてたのか. 静的多様化をつかってメタ化できるかな ?　ちょっとやってみる.

追記

やってみたらあっさり出来た. これはとっても便利だ.

#include<iostream>

template< typename DerivedMainClass >
class MainClass{
public:
    MainClass(){};
    ~MainClass(){};
    bool CallFunc(){
        return static_cast<DerivedMainClass*>(this)->CallMainFuncImple();
    }

    template< typename DerivedSubClass >
    class SubClass {
        public:
            bool CallFunc (){
                return static_cast<DerivedSubClass*>(this)->CallSubFuncImple();
            }
            SubClass(){};
            ~SubClass(){}
    };
};

class MainClassDerived : public MainClass <MainClassDerived>{
    friend class MainClass < MainClassDerived >;
    bool CallMainFuncImple(){
        std::cout << "Derived main Class" << std::endl;
    }
};

class SubClassDerived : public MainClass<MainClassDerived>::SubClass <SubClassDerived> {
    friend class SubClass < SubClassDerived >;
    bool CallSubFuncImple(){
        std::cout << "Derived sub Class" << std::endl;
    }
};

int main(){
    MainClass<MainClassDerived> *mainclass = new MainClassDerived();
    MainClass<MainClassDerived>::SubClass<SubClassDerived> *subclass = new SubClassDerived();
    mainclass->CallFunc();
    subclass->CallFunc();
return 0;}