| 1 | 
 | 
|---|
| 2 | namespace sc {
 | 
|---|
| 3 | 
 | 
|---|
| 4 | class LocalOSSContribution {
 | 
|---|
| 5 |   private:
 | 
|---|
| 6 |     double * const gmat;
 | 
|---|
| 7 |     double * const gmata;
 | 
|---|
| 8 |     double * const gmatb;
 | 
|---|
| 9 | 
 | 
|---|
| 10 |     double * const pmat;
 | 
|---|
| 11 |     double * const pmata;
 | 
|---|
| 12 |     double * const pmatb;
 | 
|---|
| 13 |   public:
 | 
|---|
| 14 |     LocalOSSContribution(double *g, double *p, double *ga, double *pa,
 | 
|---|
| 15 |                          double *gb, double *pb) :
 | 
|---|
| 16 |       gmat(g), gmata(ga), gmatb(gb), pmat(p), pmata(pa), pmatb(pb) {}
 | 
|---|
| 17 |     ~LocalOSSContribution() {}
 | 
|---|
| 18 | 
 | 
|---|
| 19 |     void set_bound(double,double) {}
 | 
|---|
| 20 | 
 | 
|---|
| 21 |     inline void cont1(int ij, int kl, double val) {
 | 
|---|
| 22 |       gmat[ij] += val*pmat[kl];
 | 
|---|
| 23 |       gmat[kl] += val*pmat[ij];
 | 
|---|
| 24 |     }
 | 
|---|
| 25 |     
 | 
|---|
| 26 |     inline void cont2(int ij, int kl, double val) {
 | 
|---|
| 27 |       val *= 0.25;
 | 
|---|
| 28 |       gmat[ij] -= val*pmat[kl];
 | 
|---|
| 29 |       gmat[kl] -= val*pmat[ij];
 | 
|---|
| 30 | 
 | 
|---|
| 31 |       gmata[ij] += val*pmata[kl];
 | 
|---|
| 32 |       gmata[kl] += val*pmata[ij];
 | 
|---|
| 33 | 
 | 
|---|
| 34 |       gmatb[ij] += val*pmatb[kl];
 | 
|---|
| 35 |       gmatb[kl] += val*pmatb[ij];
 | 
|---|
| 36 | 
 | 
|---|
| 37 |       val *= -3.0;
 | 
|---|
| 38 |       gmatb[ij] += val*pmata[kl];
 | 
|---|
| 39 |       gmatb[kl] += val*pmata[ij];
 | 
|---|
| 40 | 
 | 
|---|
| 41 |       gmata[ij] += val*pmatb[kl];
 | 
|---|
| 42 |       gmata[kl] += val*pmatb[ij];
 | 
|---|
| 43 |     }
 | 
|---|
| 44 |     
 | 
|---|
| 45 |     inline void cont3(int ij, int kl, double val) {
 | 
|---|
| 46 |       val *= 0.5;
 | 
|---|
| 47 |       gmat[ij] -= val*pmat[kl];
 | 
|---|
| 48 |       gmat[kl] -= val*pmat[ij];
 | 
|---|
| 49 | 
 | 
|---|
| 50 |       gmata[ij] += val*pmata[kl];
 | 
|---|
| 51 |       gmata[kl] += val*pmata[ij];
 | 
|---|
| 52 | 
 | 
|---|
| 53 |       gmatb[ij] += val*pmatb[kl];
 | 
|---|
| 54 |       gmatb[kl] += val*pmatb[ij];
 | 
|---|
| 55 | 
 | 
|---|
| 56 |       val *= -3.0;
 | 
|---|
| 57 |       gmata[ij] += val*pmatb[kl];
 | 
|---|
| 58 |       gmata[kl] += val*pmatb[ij];
 | 
|---|
| 59 | 
 | 
|---|
| 60 |       gmatb[ij] += val*pmata[kl];
 | 
|---|
| 61 |       gmatb[kl] += val*pmata[ij];
 | 
|---|
| 62 |     }
 | 
|---|
| 63 |     
 | 
|---|
| 64 |     inline void cont4(int ij, int kl, double val) {
 | 
|---|
| 65 |       gmat[ij] += 0.75*val*pmat[kl];
 | 
|---|
| 66 |       gmat[kl] += 0.75*val*pmat[ij];
 | 
|---|
| 67 | 
 | 
|---|
| 68 |       gmata[ij] += 0.25*val*pmata[kl];
 | 
|---|
| 69 |       gmata[kl] += 0.25*val*pmata[ij];
 | 
|---|
| 70 | 
 | 
|---|
| 71 |       gmatb[ij] += 0.25*val*pmatb[kl];
 | 
|---|
| 72 |       gmatb[kl] += 0.25*val*pmatb[ij];
 | 
|---|
| 73 | 
 | 
|---|
| 74 |       gmata[ij] -= 0.75*val*pmatb[kl];
 | 
|---|
| 75 |       gmata[kl] -= 0.75*val*pmatb[ij];
 | 
|---|
| 76 | 
 | 
|---|
| 77 |       gmatb[ij] -= 0.75*val*pmata[kl];
 | 
|---|
| 78 |       gmatb[kl] -= 0.75*val*pmata[ij];
 | 
|---|
| 79 |     }
 | 
|---|
| 80 |     
 | 
|---|
| 81 |     inline void cont5(int ij, int kl, double val) {
 | 
|---|
| 82 |       val *= 0.5;
 | 
|---|
| 83 |       gmat[ij] += val*pmat[kl];
 | 
|---|
| 84 |       gmat[kl] += val*pmat[ij];
 | 
|---|
| 85 | 
 | 
|---|
| 86 |       gmata[ij] += val*pmata[kl];
 | 
|---|
| 87 |       gmata[kl] += val*pmata[ij];
 | 
|---|
| 88 | 
 | 
|---|
| 89 |       gmatb[ij] += val*pmatb[kl];
 | 
|---|
| 90 |       gmatb[kl] += val*pmatb[ij];
 | 
|---|
| 91 | 
 | 
|---|
| 92 |       val *= -3.0;
 | 
|---|
| 93 |       gmata[ij] += val*pmatb[kl];
 | 
|---|
| 94 |       gmata[kl] += val*pmatb[ij];
 | 
|---|
| 95 | 
 | 
|---|
| 96 |       gmatb[ij] += val*pmata[kl];
 | 
|---|
| 97 |       gmatb[kl] += val*pmata[ij];
 | 
|---|
| 98 |     }
 | 
|---|
| 99 | };
 | 
|---|
| 100 | 
 | 
|---|
| 101 | class LocalOSSEnergyContribution {
 | 
|---|
| 102 |   private:
 | 
|---|
| 103 |     double * const pmat;
 | 
|---|
| 104 |     double * const pmata;
 | 
|---|
| 105 |     double * const pmatb;
 | 
|---|
| 106 | 
 | 
|---|
| 107 |   public:
 | 
|---|
| 108 |     double ec;
 | 
|---|
| 109 |     double ex;
 | 
|---|
| 110 | 
 | 
|---|
| 111 |     void set_bound(double,double) {};
 | 
|---|
| 112 |     
 | 
|---|
| 113 |     LocalOSSEnergyContribution(double *p, double *pa, double *pb) :
 | 
|---|
| 114 |       pmat(p), pmata(pa), pmatb(pb) {
 | 
|---|
| 115 |       ec=ex=0;
 | 
|---|
| 116 |     }
 | 
|---|
| 117 |     ~LocalOSSEnergyContribution() {}
 | 
|---|
| 118 | 
 | 
|---|
| 119 |     inline void cont1(int ij, int kl, double val) {
 | 
|---|
| 120 |       ec += val*pmat[ij]*pmat[kl];
 | 
|---|
| 121 |     }
 | 
|---|
| 122 |     
 | 
|---|
| 123 |     inline void cont2(int ij, int kl, double val) {
 | 
|---|
| 124 |       ex -= 0.25*val*(pmat[ij]*pmat[kl] + pmata[ij]*pmata[kl] +
 | 
|---|
| 125 |                       pmatb[ij]*pmatb[kl]);
 | 
|---|
| 126 |       ex += 0.75*val*(pmata[ij]*pmatb[kl] + pmatb[ij]*pmata[kl]);
 | 
|---|
| 127 |     }
 | 
|---|
| 128 |     
 | 
|---|
| 129 |     inline void cont3(int ij, int kl, double val) {
 | 
|---|
| 130 |       ex -= 0.5*val*(pmat[ij]*pmat[kl] + pmata[ij]*pmata[kl] +
 | 
|---|
| 131 |                       pmatb[ij]*pmatb[kl]);
 | 
|---|
| 132 |       ex += 1.5*val*(pmata[ij]*pmatb[kl] + pmatb[ij]*pmata[kl]);
 | 
|---|
| 133 |     }
 | 
|---|
| 134 |     
 | 
|---|
| 135 |     inline void cont4(int ij, int kl, double val) {
 | 
|---|
| 136 |       ec += val*pmat[ij]*pmat[kl];
 | 
|---|
| 137 |       ex -= 0.25*val*(pmat[ij]*pmat[kl] + pmata[ij]*pmata[kl] +
 | 
|---|
| 138 |                       pmatb[ij]*pmatb[kl]);
 | 
|---|
| 139 |       ex += 0.75*val*(pmata[ij]*pmatb[kl] + pmatb[ij]*pmata[kl]);
 | 
|---|
| 140 |     }
 | 
|---|
| 141 |     
 | 
|---|
| 142 |     inline void cont5(int ij, int kl, double val) {
 | 
|---|
| 143 |       ec += val*pmat[ij]*pmat[kl];
 | 
|---|
| 144 |       ex -= 0.5*val*(pmat[ij]*pmat[kl] + pmata[ij]*pmata[kl] +
 | 
|---|
| 145 |                       pmatb[ij]*pmatb[kl]);
 | 
|---|
| 146 |       ex += 1.5*val*(pmata[ij]*pmatb[kl] + pmatb[ij]*pmata[kl]);
 | 
|---|
| 147 |     }
 | 
|---|
| 148 | };
 | 
|---|
| 149 | 
 | 
|---|
| 150 | class LocalOSSGradContribution {
 | 
|---|
| 151 |   private:
 | 
|---|
| 152 |     double * const pmat;
 | 
|---|
| 153 |     double * const pmata;
 | 
|---|
| 154 |     double * const pmatb;
 | 
|---|
| 155 | 
 | 
|---|
| 156 |   public:
 | 
|---|
| 157 |     LocalOSSGradContribution(double *p, double *pa, double *pb) :
 | 
|---|
| 158 |       pmat(p), pmata(pa), pmatb(pb) {}
 | 
|---|
| 159 |     ~LocalOSSGradContribution() {}
 | 
|---|
| 160 | 
 | 
|---|
| 161 |     inline double cont1(int ij, int kl) {
 | 
|---|
| 162 |       return pmat[ij]*pmat[kl] +
 | 
|---|
| 163 |         0.5*(pmata[ij]*pmat[kl] + pmat[ij]*pmata[kl] +
 | 
|---|
| 164 |              pmatb[ij]*pmat[kl] + pmat[ij]*pmatb[kl]) +
 | 
|---|
| 165 |         0.25*(pmata[ij]*pmata[kl] + pmatb[ij]*pmatb[kl] +
 | 
|---|
| 166 |               pmata[ij]*pmatb[kl] + pmatb[ij]*pmata[kl]);
 | 
|---|
| 167 |     }
 | 
|---|
| 168 | 
 | 
|---|
| 169 |     inline double cont2(int ij, int kl) {
 | 
|---|
| 170 |       return pmat[ij]*pmat[kl] +
 | 
|---|
| 171 |         0.5*(pmata[ij]*pmat[kl] + pmat[ij]*pmata[kl] +
 | 
|---|
| 172 |              pmatb[ij]*pmat[kl] + pmat[ij]*pmatb[kl] +
 | 
|---|
| 173 |              pmata[ij]*pmata[kl] + pmatb[ij]*pmatb[kl] -
 | 
|---|
| 174 |              pmata[ij]*pmatb[kl] - pmatb[ij]*pmata[kl]);
 | 
|---|
| 175 |     }
 | 
|---|
| 176 | };
 | 
|---|
| 177 | 
 | 
|---|
| 178 | }
 | 
|---|
| 179 | 
 | 
|---|