~ubuntu-branches/ubuntu/trusty/yade/trusty

« back to all changes in this revision

Viewing changes to pkg/dem/Ip2_ElastMat.cpp

  • Committer: Package Import Robot
  • Author(s): Anton Gladky, cf3f8d9
  • Date: 2013-10-30 20:56:33 UTC
  • mfrom: (20.1.9 sid)
  • Revision ID: package-import@ubuntu.com-20131030205633-1f01r7hjce17d723
Tags: 1.05.0-2
[cf3f8d9] Pass -ftrack-macro-expansion=0 only if gcc>=4.8. (Closes: #726009)

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#include "Ip2_ElastMat.hpp"
 
2
 
 
3
#include <yade/pkg/common/NormShearPhys.hpp>
 
4
#include <yade/pkg/common/NormShearPhys.hpp>
 
5
#include <yade/pkg/dem/DemXDofGeom.hpp>
 
6
 
 
7
YADE_PLUGIN((Ip2_ElastMat_ElastMat_NormPhys)(Ip2_ElastMat_ElastMat_NormShearPhys));
 
8
 
 
9
 
 
10
void Ip2_ElastMat_ElastMat_NormPhys::go( const shared_ptr<Material>& b1
 
11
                                        , const shared_ptr<Material>& b2
 
12
                                        , const shared_ptr<Interaction>& interaction)
 
13
{
 
14
        if(interaction->phys) return;
 
15
        const shared_ptr<ElastMat>& mat1 = YADE_PTR_CAST<ElastMat>(b1);
 
16
        const shared_ptr<ElastMat>& mat2 = YADE_PTR_CAST<ElastMat>(b2);
 
17
        Real Ea         = mat1->young;
 
18
        Real Eb         = mat2->young;
 
19
        interaction->phys = shared_ptr<NormPhys>(new NormPhys());
 
20
        const shared_ptr<NormPhys>& phys = YADE_PTR_CAST<NormPhys>(interaction->phys);
 
21
        Real Kn;
 
22
        const GenericSpheresContact* geom=dynamic_cast<GenericSpheresContact*>(interaction->geom.get());
 
23
        if (geom) {
 
24
                Real Ra,Rb;//Vector3r normal;
 
25
                Ra=geom->refR1>0?geom->refR1:geom->refR2;
 
26
                Rb=geom->refR2>0?geom->refR2:geom->refR1;
 
27
                //harmonic average of the two stiffnesses when (Ri.Ei/2) is the stiffness of a contact point on sphere "i"
 
28
                Kn = 2*Ea*Ra*Eb*Rb/(Ea*Ra+Eb*Rb);
 
29
        } else {
 
30
                Kn = 2*Ea*Eb/(Ea+Eb);
 
31
        }
 
32
        phys->kn = Kn;
 
33
};
 
34
 
 
35
 
 
36
void Ip2_ElastMat_ElastMat_NormShearPhys::go( const shared_ptr<Material>& b1
 
37
                                        , const shared_ptr<Material>& b2
 
38
                                        , const shared_ptr<Interaction>& interaction)
 
39
{
 
40
        if(interaction->phys) return;
 
41
        const shared_ptr<ElastMat>& mat1 = YADE_PTR_CAST<ElastMat>(b1);
 
42
        const shared_ptr<ElastMat>& mat2 = YADE_PTR_CAST<ElastMat>(b2);
 
43
        Real Ea         = mat1->young;
 
44
        Real Eb         = mat2->young;
 
45
        Real Va         = mat1->poisson;
 
46
        Real Vb         = mat2->poisson;
 
47
        interaction->phys = shared_ptr<NormShearPhys>(new NormShearPhys());
 
48
        const shared_ptr<NormShearPhys>& phys = YADE_PTR_CAST<NormShearPhys>(interaction->phys);
 
49
        Real Kn=0.0, Ks=0.0;
 
50
        GenericSpheresContact* geom=dynamic_cast<GenericSpheresContact*>(interaction->geom.get());
 
51
        if (geom) {
 
52
                Real Ra,Rb;//Vector3r normal;
 
53
                Ra=geom->refR1>0?geom->refR1:geom->refR2;
 
54
                Rb=geom->refR2>0?geom->refR2:geom->refR1;
 
55
                //harmonic average of the two stiffnesses when (Ri.Ei/2) is the stiffness of a contact point on sphere "i"
 
56
                Kn = 2*Ea*Ra*Eb*Rb/(Ea*Ra+Eb*Rb);
 
57
                Ks = 2*Ea*Ra*Va*Eb*Rb*Vb/(Ea*Ra*Va+Eb*Rb*Vb);
 
58
        } else {
 
59
                Kn = 2*Ea*Eb/(Ea+Eb);
 
60
                Kn = 2*Ea*Va*Eb*Vb/(Ea*Va+Eb*Vb);
 
61
        }
 
62
        phys->kn = Kn;
 
63
        phys->ks = Ks;
 
64
};