ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/perturbations/UniformGradient.cpp
(Generate patch)

Comparing trunk/src/perturbations/UniformGradient.cpp (file contents):
Revision 2026 by gezelter, Wed Oct 22 12:23:59 2014 UTC vs.
Revision 2047 by gezelter, Thu Dec 11 16:16:43 2014 UTC

# Line 58 | Line 58 | namespace OpenMD {
58    }
59    
60    void UniformGradient::initialize() {
61 <    if (simParams->haveUniformGradient()) {
62 <      doUniformGradient = true;
63 <      std::vector<RealType> pv = simParams->getUniformGradient();            
64 <      if (pv.size() != 5) {
61 >
62 >    bool haveA = false;
63 >    bool haveB = false;
64 >    bool haveG = false;
65 >    
66 >    if (simParams->haveUniformGradientDirection1()) {
67 >      std::vector<RealType> d1 = simParams->getUniformGradientDirection1();
68 >      if (d1.size() != 3) {
69          sprintf(painCave.errMsg,
70 <                "UniformGradient: Incorrect number of parameters specified.\n"
71 <                "\tthere should be 5 parameters, but %lu were specified.\n", pv.size());
70 >                "uniformGradientDirection1: Incorrect number of parameters\n"
71 >                "\tspecified. There should be 3 parameters, but %lu were\n"
72 >                "\tspecified.\n", d1.size());
73          painCave.isFatal = 1;
74          simError();      
75        }
76 <      pars_.a = pv[0];
77 <      pars_.b = pv[1];
78 <      pars_.c = pv[2];
74 <      pars_.alpha = pv[3];
75 <      pars_.beta = pv[4];
76 >      a_.x() = d1[0];
77 >      a_.y() = d1[1];
78 >      a_.z() = d1[2];
79  
80 <      Grad_(0,0) = pars_.alpha;
81 <      Grad_(0,1) = pars_.a;
82 <      Grad_(0,2) = pars_.b;
80 >      a_.normalize();
81 >      haveA = true;      
82 >    }
83 >    
84 >    if (simParams->haveUniformGradientDirection2()) {
85 >      std::vector<RealType> d2 = simParams->getUniformGradientDirection2();
86 >      if (d2.size() != 3) {
87 >        sprintf(painCave.errMsg,
88 >                "uniformGradientDirection2: Incorrect number of parameters\n"
89 >                "\tspecified. There should be 3 parameters, but %lu were\n"
90 >                "\tspecified.\n", d2.size());
91 >        painCave.isFatal = 1;
92 >        simError();      
93 >      }
94 >      b_.x() = d2[0];
95 >      b_.y() = d2[1];
96 >      b_.z() = d2[2];
97 >
98 >      b_.normalize();
99 >      haveB = true;      
100 >    }
101 >
102 >    if (simParams->haveUniformGradientStrength()) {
103 >      g_ = simParams->getUniformGradientStrength();
104 >      haveG = true;
105 >    }
106 >
107 >    if (haveA && haveB && haveG) {
108 >      doUniformGradient = true;
109 >      cpsi_ = dot(a_, b_);
110 >      
111 >      Grad_(0,0) = 2.0 * (a_.x()*b_.x() - cpsi_ / 3.0);
112 >      Grad_(0,1) = a_.x()*b_.y() + a_.y()*b_.x();
113 >      Grad_(0,2) = a_.x()*b_.z() + a_.z()*b_.x();
114        Grad_(1,0) = Grad_(0,1);
115 <      Grad_(1,1) = pars_.beta;
116 <      Grad_(1,2) = pars_.c;
115 >      Grad_(1,1) = 2.0 * (a_.y()*b_.y() - cpsi_ / 3.0);
116 >      Grad_(1,2) = a_.y()*b_.z() + a_.z()*b_.y();
117        Grad_(2,0) = Grad_(0,2);
118        Grad_(2,1) = Grad_(1,2);
119 <      Grad_(2,2) = - (Grad_(0,0) + Grad_(1,1));
120 <    }  
119 >      Grad_(2,2) = 2.0 * (a_.z()*b_.z() - cpsi_ / 3.0);
120 >
121 >      Grad_ *= g_ / 2.0;
122 >
123 >    } else {
124 >      if (!haveA) {
125 >        sprintf(painCave.errMsg,
126 >                "UniformGradient: uniformGradientDirection1 not specified.\n");
127 >        painCave.isFatal = 1;
128 >        simError();
129 >      }
130 >      if (!haveB) {
131 >        sprintf(painCave.errMsg,
132 >                "UniformGradient: uniformGradientDirection2 not specified.\n");
133 >        painCave.isFatal = 1;
134 >        simError();
135 >      }
136 >      if (!haveG) {
137 >        sprintf(painCave.errMsg,
138 >                "UniformGradient: uniformGradientStrength not specified.\n");
139 >        painCave.isFatal = 1;
140 >        simError();
141 >      }
142 >    }    
143 >    
144      int storageLayout_ = info_->getSnapshotManager()->getStorageLayout();
145      if (storageLayout_ & DataStorage::dslParticlePot) doParticlePot = true;
146      initialized = true;

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines