+++ /dev/null
-#include <iostream>
-#include "McCaskill.hpp"
-//#include "energy_param3.hpp"
-#include "Util.hpp"
-#include <cstring>
-
-namespace MXSCARNA {
-energy_param McCaskill::energyParam;
-
-float *McCaskill::exphairpin;
-float McCaskill::expninio[5][32];
-float McCaskill::expdangle5[6][4];
-float McCaskill::expdangle3[6][4];
-float McCaskill::expinternalLoop[31];
-float McCaskill::expbulge[31];
-char McCaskill::exptriLoop[2][6];
-float McCaskill::exptriLoopEnergy[2];
-char McCaskill::exptetraLoop[30][7];
-float McCaskill::exptetraLoopEnergy[30];
-float McCaskill::expStack[6][6];
-float McCaskill::expTstackH[6][16];
-float McCaskill::expTstackI[6][16];
-float McCaskill::expint11[6][6][16];
-float McCaskill::expint21[6][6][16][4];
-float McCaskill::expint22[6][6][16][16];
-float McCaskill::expMLclosing;
-float McCaskill::expMLintern[8];
-float McCaskill::expTermAU;
-float McCaskill::expMLbase[31];
-
-const int McCaskill::TURN = 3;
-const float McCaskill::GASCONST = 1.98717;
-const float McCaskill::T = 37 + 273.15;
-const int McCaskill::MAXLOOP = 30;
-const int McCaskill::TETRA_ENTH37 = -400;
-const int McCaskill::NBPAIRS = 7;
-const int McCaskill::SCALE = 10;
-
-
-void
-McCaskill::
-calcPartitionFunction()
-{
- initParameter();
- Inside();
- Outside();
- convertProbability();
-
-/*
- for(int i = 0; i < n_seq; i++) {
- for(int j = i; j < n_seq; j++) {
- cout << getProb(i, j) << " ";
- }
- cout << endl;
- }
-*/
-}
-
-void
-McCaskill::
-convertProbability()
-{
- float *pPointer = p.getPointer(0, 0);
- float *abPointer = ab.getPointer(0,0);
- for(int i = 0; i < n_seq; i++) {
- for(int j = i; j < n_seq; j++) {
- *pPointer += *abPointer++;
- *pPointer++ = EXP(*pPointer);
- }
- }
-}
-
-void
-McCaskill::
-Inside()
-{
-
- for (int i = n_seq - TURN - 2; i >= 0; i--) {
- float *abPointer = ab.getPointer(i, i + TURN + 1);
- float *am1Pointer = am1.getPointer(i, i + TURN + 1);
- float *amPointer = am.getPointer(i, i + TURN + 1);
- float *q1Pointer = q1.getPointer(i, i + TURN + 1);
- float *aPointer = a.getPointer(i, i + TURN + 1);
- int *typePointer = typeMat.getPointer(i, i+TURN+1);
- for (int j = i + TURN + 1; j < n_seq; j++) {
- int tmpType = *typePointer++;
- *abPointer++ = compQb(i, j, tmpType);
- am1v.ref(i,j) = *am1Pointer++ = compQm1(i,j, tmpType);
- amv.ref(i,j) = *amPointer++ = compQm(i,j);
- q1v.ref(i,j) = *q1Pointer++ = compQ1(i,j, tmpType);
- *aPointer++ = compQ(i,j);
- }
- }
-}
-
-inline float
-McCaskill::
-compQb(int i, int j, int tmpType)
-{
-
- float tmpAb;
- int type = tmpType;
- int u = j - i - 1;
-
- if(Beta::isCanonicalReducedPairCode(type)) {
- // hairpin energy
- assert(u >= 3);
- tmpAb = expHairpinEnergy(type, u, i + 1, j - 1);
-
- // base pair, bulge, interior loop energy
- for(int h = i + 1; h <= MIN(i + MAXLOOP + 1, j - TURN - 2); h++) {
- int u1 = h-i-1;
- int max = MAX(h + TURN + 1, j-1-MAXLOOP+u1);
- float *abPointer = ab.getPointer(h, max - 1);
- const int *typeMatPointer = typeMat.getPointer(h, max);
-
- for(int l = max; l < j; l++) {
- int type2 = *typeMatPointer++;
- abPointer++;
- if(!Beta::isCanonicalReducedPairCode(type2)) continue;
-
- assert(h >= 0 && h < n_seq && l >= 0 && l < n_seq);
- type2 = Beta::flipReducedPairCode(type2);
- assert(h-i-1 >= 0); assert(j-l-1 >= 0);
- float loopE = *abPointer;
- loopE += expLoopEnergy(u1, j-l-1, tmpType, type2, i+1, j-1, h-1, l+1);
- tmpAb = LOG_ADD(tmpAb, loopE);
- }
- }
-
- // multi loop
- float tmpQm = IMPOSSIBLE;
- float *amPointer = am.getPointer(i + 1, j - TURN - 3);
- float *am1Pointer = am1v.getPointer(j-TURN-2, j-1);
- for(int h = j - TURN - 2; h >= i + TURN + 3; h--) {
- assert(h >= 0 && h < n_seq);
- float tmpScore = *amPointer--;
- tmpScore += *am1Pointer--;
- tmpQm = LOG_ADD(tmpQm, tmpScore);
- }
- tmpQm += expMLclosing + expMLintern[type];
- tmpQm += endStemScore(i, j);
- tmpAb = LOG_ADD(tmpAb, tmpQm);
- }
- else {
- tmpAb = IMPOSSIBLE;
- }
- return tmpAb;
-}
-
-//F = a:ML_closing + b:ML_intern*k + c:ML_BASE*u
-
-inline float
-McCaskill::
-compQm1(int i, int j, int tmpType)
-{
- float tmpQm1 = IMPOSSIBLE;
-
- int l = j;
- if (i + TURN + 1 <= l) {
- int type = typeMat.ref(i,l);
- if(Beta::isCanonicalReducedPairCode(type)) {
- float tmpScore = ab.ref(i,l);
- tmpScore += beginStemScore(i, l);
- //tmpScore += expMLbase[1]*(j-l) + expMLintern[type];
- tmpScore += expMLintern[tmpType];
- tmpQm1 = LOG_ADD(tmpQm1, tmpScore);
- }
- }
- if(i < j) {
- tmpQm1 = LOG_ADD(tmpQm1, am1.ref(i,j-1));
- }
-
- return tmpQm1;
-}
-
-inline float
-McCaskill::
-compQm(int i, int j)
-{
- float tmpQm = IMPOSSIBLE;
- float *amPointer = am.getPointer(i,j-TURN-2);
- float *am1Pointer = am1v.getPointer(j-TURN-1, j);
- for(int h = j - TURN - 1; h >= i ; h--) {
- float tmpScore = 0;
- float tmpAm1 = *am1Pointer--;
-
- tmpScore += tmpAm1;
- tmpQm = LOG_ADD(tmpQm, tmpScore);
- tmpScore = *amPointer--;
- tmpScore += tmpAm1;
- tmpQm = LOG_ADD(tmpQm, tmpScore);
- }
-
- return tmpQm;
-}
-
-inline float
-McCaskill::
-compQ1(int i, int j, int tmpType)
-{
- float tmpQ1 = IMPOSSIBLE;
-
- if(Beta::isCanonicalReducedPairCode(tmpType)) {
- float tmpScore = ab.ref(i, j);
- tmpScore += beginStemScore(i, j);
- tmpQ1 = LOG_ADD(tmpQ1, tmpScore);
- }
- tmpQ1 = LOG_ADD(tmpQ1, q1.ref(i, j - 1));
-
- return tmpQ1;
-}
-
-inline float
-McCaskill::
-compQ(int i, int j)
-{
- float tmpQ = 0;
- tmpQ = LOG_ADD(tmpQ, q1.ref(i,j));
-
- float *aPointer = a.getPointer(i,j-TURN-2);
- float *q1Pointer = q1v.getPointer(j-TURN-1, j);
- for(int h = j - TURN - 1; h >= i + 1; h--) {
- float tmpScore = *aPointer--;
- tmpScore += *q1Pointer--;
- tmpQ = LOG_ADD(tmpQ, tmpScore);
- }
-
- return tmpQ;
-}
-
-inline float
-McCaskill::
-beginStemScore(const int i, const int j) const
-{
- float temp = 0;
- int type = typeMat.ref(i,j);
-
- if(0 < i) { temp += expdangle5[type][numSeq[i-1]]; }
- if(j < n_seq-1) { temp += expdangle3[type][numSeq[j+1]]; }
- if(type != Beta::REDUCED_CG_CODE && type != Beta::REDUCED_GC_CODE) { temp += expTermAU; }
- return temp;
-}
-
-inline float
-McCaskill::
-endStemScore(const int i, const int j) const
-{
- float temp = 0;
- int type = typeMat.ref(i,j);
-
- type = Beta::flipReducedPairCode(type);
-
- if(i < n_seq-1) { temp += expdangle3[type][numSeq[i+1]]; }
- if(j > 0) { temp += expdangle5[type][numSeq[j-1]]; }
- if(type != Beta::REDUCED_CG_CODE && type != Beta::REDUCED_GC_CODE) { temp += expTermAU; }
- return temp;
-}
-
-inline float
-McCaskill::
-compP(int h, int l, int tmpType)
-{
- float prob = IMPOSSIBLE;
-
- int type = tmpType;
- if(Beta::isCanonicalReducedPairCode(type)) {
-
- /* exterior loop */
- float tmp_p = 0;
- tmp_p -= a.ref(0,n_seq-1);
- if(0 < h) {
- tmp_p += a.ref(0,h-1);
- }
- if(l < n_seq-1) {
- tmp_p += a.ref(l+1, n_seq-1);
- }
- tmp_p += beginStemScore(h, l);
- prob = LOG_ADD(prob, tmp_p);
-
- assert(IMPOSSIBLE <= prob && prob <= 0);
-
- /* internal loop */
- tmp_p = IMPOSSIBLE;
- int tt = Beta::flipReducedPairCode(tmpType);
- int max = MAX(0,h-MAXLOOP-1);
- for(int i = max; i <= h - 1; i++) {
- float min = MIN(l+MAXLOOP-h+i+2, n_seq-1);
- int *typeMatPointer = typeMat.getPointer(i,l+1);
- float *pPointer = p.getPointer(i,l);
- for(int j = l + 1; j <= min; j++) {
- int type2 = *typeMatPointer++;
- pPointer++;
- if(!Beta::isCanonicalReducedPairCode(type2)) continue;
- assert(i >= 0 && i < n_seq && j >= 0 && j < n_seq);
-
- float tmpScore = *pPointer;
- tmpScore += expLoopEnergy(h-i-1, j-l-1, type2, tt, i+1, j-1, h-1, l+1);
- tmp_p = LOG_ADD(tmp_p, tmpScore);
- }
- }
- prob = LOG_ADD(prob, tmp_p);
- assert(IMPOSSIBLE <= prob && prob <= 0);
-
- /* multi loop */
- tmp_p = IMPOSSIBLE;
- float tmp_begin = beginStemScore(h, l);
- float *q1Pointer = q1v.getPointer(0, l);
- float *am1Pointer = am1v.getPointer(0, l);
- float *amPointer = amv.getPointer(1,h-1);
- for(int i = 0; i <= h-TURN-1; i++) {
- float tmpq1 = *q1Pointer++;
- float tmpam = *amPointer++;
- float tmpScore = *am1Pointer++;
-
- tmpScore += tmpam;
- tmpScore += tmp_begin;
- tmpScore += expMLclosing + expMLintern[tt];
- tmp_p = LOG_ADD(tmp_p, tmpScore);
-
- tmpScore = tmpq1;
- tmpScore += tmpam;
- tmpScore += tmp_begin;
- tmpScore += expMLclosing + expMLintern[tt];
- tmp_p = LOG_ADD(tmp_p, tmpScore);
-
- tmpScore = tmpq1;
- tmpScore += tmp_begin;
- tmpScore += expMLclosing + expMLintern[tt];
- tmp_p = LOG_ADD(tmp_p, tmpScore);
- }
-
- assert(IMPOSSIBLE <= tmp_p && tmp_p <= 0);
- prob = LOG_ADD(prob, tmp_p);
-
- tmp_p = IMPOSSIBLE;
- for(int i = h-TURN; i <= h-1; i++) {
- if(i >= 0) {
- float tmpScore = q1.ref(i,l);
- tmpScore += tmp_begin;
- tmpScore += expMLclosing + expMLintern[tt];
- tmp_p = LOG_ADD(tmp_p, tmpScore);
- }
- }
- assert(IMPOSSIBLE <= tmp_p && tmp_p <= 0);
- prob = LOG_ADD(prob, tmp_p);
- }
- else {
- prob = IMPOSSIBLE;
- }
-
- return prob;
-}
-
-inline float
-McCaskill::
-compPm(int i, int l)
-{
- float tmpPm = IMPOSSIBLE;
-
- int *typeMatPointer = typeMat.getPointer(i,n_seq-1);
- float *pPointer = p.getPointer(i,n_seq);
- float *amPointer = am.getPointer(l+1,n_seq-1);
- float *abPointer = ab.getPointer(i, n_seq);
- for(int j = n_seq - 1; j >= l + TURN + 1; j--) {
- int type = *typeMatPointer--;
- pPointer--;
- amPointer--;
- abPointer--;
- if(Beta::isCanonicalReducedPairCode(type)) {
- float tmp = *pPointer;
- tmp += *amPointer;
- tmp += endStemScore(i, j);
- tmpPm = LOG_ADD(tmpPm, tmp);
- }
- }
- tmpPm += expMLintern[1];
-
- return tmpPm;
-}
-
-inline float
-McCaskill::
-compPm1(int i, int l)
-{
- float tmpPm1 = IMPOSSIBLE;
-
- int j = l + 1;
- if(j <= n_seq-1) {
- int type = typeMat.ref(i,j);
- if(Beta::isCanonicalReducedPairCode(type)) {
- float tmp = p.ref(i,j);
- tmp += endStemScore(i, j);
- tmpPm1 = tmp;
- }
- tmpPm1 += expMLintern[1];
- }
- if(l+1 <= n_seq - 1) {
- tmpPm1 = LOG_ADD(tmpPm1, am1.ref(i, l+1));
- }
-
- return tmpPm1;
-}
-
-void
-McCaskill::
-Outside()
-{
- for(int h = 0; h <= n_seq - TURN - 2; h++) {
- float *pPointer = p.getPointer(h, n_seq-1);
- float *q1Pointer = q1.getPointer(h, n_seq-1);
- float *am1Pointer = am1.getPointer(h, n_seq-1);
- int *typePointer = typeMat.getPointer(h, n_seq-1);
- for(int l = n_seq-1; l >= h + TURN + 1; l--) {
- int tmpType = *typePointer--;
- pv.ref(h,l) = *pPointer-- = compP(h,l,tmpType);
- q1v.ref(h,l) = *q1Pointer-- = compPm(h,l);
- am1v.ref(h,l) = *am1Pointer-- = compPm1(h,l);
-
- assert(p.ref(h,l) <= 0);
- }
- }
-}
-
-void
-McCaskill::
-printProbMat()
-{
- int m = 0;
- for(int i = 0; i < n_seq; i++) cout << " " << seq[i];
- cout << endl;
- for(int i = 0; i < n_seq; i++) {
- if(m < n_seq) {
- cout << seq[m];
- }
- for(int j = 0; j <= i-1; j++) {
- if(j != i-1) cout << " ";
- else cout << " ";
- }
- if(i != 0 && i != n_seq-1) {
- cout << "\\";
- }
-
- for(int j = i; j < n_seq; j++) {
- if(p.ref(i,j) > 0.01) {
-
- int type = Beta::getReducedPairCode(numSeq[i], numSeq[j]);
-
- if(!Beta::isCanonicalReducedPairCode(type)) {
- cout << "\n" << seq[i] << " " << seq[j] << " " << exp(p.ref(i,j)) << endl;
- }
-
- if(j != n_seq-1) {
- cout << "* ";
- }
- else {
- cout << "*";
- }
-
- }
- else {
-
- if(j != n_seq-1) {
- cout << " ";
- }
- else {
- cout << " ";
- }
-
- }
-
- }
- if(m < n_seq) {
- cout << seq[m++] << endl;
- }
- if(i == n_seq - 1) cout << endl;
- }
- for(int i = 0; i < n_seq; i++) cout << " " << seq[i];
- cout << endl;
-}
-
-void
-McCaskill::
-initParameter()
-{
- float GT;
- float RT_KCAL_MOL = McCaskill::T*McCaskill::GASCONST;
- int len = 31;
-
- for(int i = 0; i < len; i++) {
- GT = energyParam.getHairpin(i);
- exphairpin[i] = -GT*10/RT_KCAL_MOL;
- }
-
- for (int i = 0; i < len; i++) {
- GT = energyParam.getBulge(i);
- expbulge[i] = -GT*10/RT_KCAL_MOL;
- GT = energyParam.getInternalLoop(i);
- expinternalLoop[i] = -GT*10/RT_KCAL_MOL;
- }
- expinternalLoop[2] = -80*10/RT_KCAL_MOL; /* special case of size 2 interior loops (single mismatch) */
-
- // float lxc = energy_param3::lxc37;
- for (int i = 31; i < n_seq; i++) {
- GT = energyParam.getHairpin(30) + (107.856*LOG((float)i/30));
- exphairpin[i] = -GT*10/RT_KCAL_MOL;
- }
-
- for(int i = 0; i < 5; i++) {
- GT = energyParam.getNinio(i);
- for(int j = 0; j <= MAXLOOP; j++) {
- expninio[i][j] = -MIN(energyParam.getMaxNinio(), j*GT)*10/RT_KCAL_MOL;
- }
- }
-
- for(int i = 0; i < 30; i++) {
- GT = energyParam.getTetraLoopEnergy(i);
- exptetraLoopEnergy[i] = -GT*10/RT_KCAL_MOL;
- }
-
- /*no parameters for Triloop*/
- for(int i = 0; i < 2; i++) {
- GT = 0;
- exptriLoopEnergy[i] = -GT*10/RT_KCAL_MOL;
- }
-
- GT = energyParam.getMLclosing();
- expMLclosing = -GT*10/RT_KCAL_MOL;
-
- for(int i = 0; i <= NBPAIRS; i++) {
- GT = energyParam.getMLintern();
- expMLintern[i] = -GT*10/RT_KCAL_MOL;
- }
-
- expTermAU = -energyParam.getTerminalAU()*10/RT_KCAL_MOL;
-
- GT = energyParam.getMLBase();
- for(int i = 0; i < len; i++) {
- expMLbase[i] = -GT*10*(float)i/RT_KCAL_MOL;
- }
-
- /*
- if danlges = 0 just set their energy to 0,
- don't let dangle energyies become > 0 (at large temps)
- but make sure go smoothly to 0
- */
- for(int i = 0; i < 6; i++) {
- for(int j =0; j < 4; j++) {
- GT = energyParam.getDangle5(i,j);
- expdangle5[i][j] = -GT*10/RT_KCAL_MOL;
- GT = energyParam.getDangle3(i,j);
- expdangle3[i][j] = -GT*10/RT_KCAL_MOL;
- }
- }
-
- /* stacking energies */
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 6; j++) {
- GT = energyParam.getStack(i,j);
- expStack[i][j] = -GT*10/RT_KCAL_MOL;
- }
- }
-
- /* mismatch energies */
- for (int i = 0; i < 6; i++) {
- for (int j = 0; j < 16; j++) {
- GT = energyParam.getTstackI(i, j);
- // cout << i << " " << " " << j << " " << GT << endl;
- expTstackI[i][j] = -GT*10/RT_KCAL_MOL;
- GT = energyParam.getTstackH(i, j);
- expTstackH[i][j] = -GT*10/RT_KCAL_MOL;
- }
- }
-
- /* interior loops of length 2*/
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 6; j++) {
- for(int k = 0; k < 16; k++) {
- GT = energyParam.getInt11(i, j, k);
- expint11[i][j][k] = -GT*10/RT_KCAL_MOL;
- }
- }
- }
-
- /* interior 2*1 loops */
- for(int i = 0; i < 6; i++) {
- for(int j =0; j < 6; j++) {
- for(int k = 0; k < 16; k++) {
- for(int l = 0; l < 4; l++) {
- GT = energyParam.getInt21(i,j,k,l);
- expint21[i][j][k][l] = -GT*10/RT_KCAL_MOL;
- }
- }
- }
- }
-
- /* interior 2*2 loops */
- for (int i = 0; i < 6; i++) {
- for(int j = 0; j < 6; j++) {
- for(int k = 0; k < 16; k++) {
- for(int l = 0; l < 16; l++) {
- GT = energyParam.getInt22(i,j,k,l);
- expint22[i][j][k][l] = -GT*10/RT_KCAL_MOL;
- }
- }
- }
- }
-}
-
-
-inline float
-McCaskill::
-expHairpinEnergy(const int type, const int l, const int i, const int j)
-{
- float q;
- int k;
-
-// assert(l >= 0);
- q = exphairpin[l];
-
- if(l == 4) {
- char temp_seq[7];
-
- for(int iter = i - 1; iter < i + 5; iter++) {
- temp_seq[iter - (i-1)] = seq[iter];
- }
- temp_seq[6] = '\0';
-
- for(k = 0; k < 30; k++) {
- if(strcmp(temp_seq, energyParam.getTetraLoop(k)) == 0) break;
- }
- if(k != 30) {
- q += exptetraLoopEnergy[k];
- }
- }
- if(l == 3) {
-
- /* no triloop bonus
- char temp_seq[6];
-
- for(int iter = i - 1; iter < i + 4; iter++) {
- temp_seq[iter - (i-1)] = seq[iter];
- }
- temp_seq[6] = '\0';
- for(k = 0; k < 2; k++) {
- if(strcmp(temp_seq, energyParam.getTriLoop(k)) == 0) break;
- }
- if(k != 2) {
- q *= exptriLoopEnergy[k];
- }
- */
-
- if(type != Beta::REDUCED_CG_CODE && type != Beta::REDUCED_GC_CODE) q += expTermAU;
- }
- else {
- int type2 = Beta::getPairCode(numSeq[i], numSeq[j]);
- q += expTstackH[type][type2];
- }
-
- return q;
-}
-
-
-inline float
-McCaskill::
-expLoopEnergy(int u1, int u2, int type, int type2,
- int si1, int sj1, int sp1, int sq1)
-{
- float z = 0;
-
- if((u1 == 0) && (u2 == 0)) { z = expStack[type][type2]; }
- else {
- if((u1 == 0) || (u2 == 0)) {
- int u;
- if(u1 == 0) { u = u2; }
- else { u = u1; }
- z = expbulge[u];
- if(u1 + u2 == 1) z += expStack[type][type2];
- else {
- if (type != Beta::REDUCED_CG_CODE && type != Beta::REDUCED_GC_CODE) z += expTermAU;
- if (type2 != Beta::REDUCED_CG_CODE && type2 != Beta::REDUCED_GC_CODE) z += expTermAU;
- }
- }
- else {
- if(u1 + u2 == 2) {
- z = expint11[type][type2][Beta::getPairCode(numSeq[si1], numSeq[sj1])];
- }
- else if((u1 == 1) && (u2 == 2)) {
- z = expint21[type][type2][Beta::getPairCode(numSeq[si1], numSeq[sj1])][numSeq[sq1]];
- }
- else if((u1 == 2) && (u2 == 1)) {
- z = expint21[type2][type][Beta::getPairCode(numSeq[sq1], numSeq[sp1])][numSeq[si1]];
- }
- else if((u1 == 2) && (u2 == 2)) {
- z = expint22[type][type2][Beta::getPairCode(numSeq[si1], numSeq[sj1])][Beta::getPairCode(numSeq[sp1], numSeq[sq1])];
- }
- else {
- z = expinternalLoop[u1 + u2] +
- expTstackI[type][Beta::getPairCode(numSeq[si1], numSeq[sj1])]
- + expTstackI[type2][Beta::getPairCode(numSeq[sq1], numSeq[sp1])];
- z += expninio[2][abs(u1-u2)];
- }
- }
- }
-
- return z;
-}
-
-void
-McCaskill::
-printExpEnergy()
-{
- cout << "exphairpin:" << endl;
- for(int i = 0; i < 31; i++) {
- cout << exphairpin[i] << endl;
- }
-
- cout << "expninio[5][32]:" << endl;
- for(int i = 0; i < 5; i++) {
- for(int j = 0; j < 32; j++) {
- cout << expninio[i][j] << " ";
- }
- cout << endl;
- }
-
- cout << "expdangle5[6][4]:" << endl;
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 4; j++) {
- cout << expdangle5[i][j] << " ";
- }
- cout << endl;
- }
-
- cout << "expdangle3[6][4]:" << endl;
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 4; j++) {
- cout << expdangle3[i][j] << " ";
- }
- cout << endl;
- }
-
- cout << "expinternalLoop[31]:" << endl;
- for(int i = 0; i < 31; i++) {
- cout << i << ":" << expinternalLoop[i] << endl;
- }
- cout << "expbulge[31]:" << endl;
- for(int i = 0; i < 31; i++) {
- cout << i << ":" << expbulge[i] << endl;
- }
-
- cout << "exptriLoopEnergy[2]:" << endl;
- for(int i = 0; i < 2; i++) {
- cout << i << ":" << exptriLoopEnergy[i] << endl;
- }
-
- cout << "exptetraLoopEnergy[15]" << endl;
- for(int i = 0; i < 15; i++) {
- cout << i << ":" << exptetraLoopEnergy[i] << endl;
- }
-
- cout << "expStack[6][6]:" << endl;
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 6; j++) {
- cout << expStack[i][j] << " ";
- }
- cout << endl;
- }
-
- cout << "expTstackH[6][16]:" << endl;
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 16; j++) {
- cout << expTstackH[i][j] << " ";
- }
- cout << endl;
- }
-
- cout << "expTstackI[6][16]:" << endl;
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 16; j++) {
- cout << expTstackI[i][j] << " ";
- }
- cout << endl;
- }
-
- cout << "expMLclosing=" << expMLclosing << endl;
- cout << "expMLintern:" << endl;
- for(int i = 0; i < 8; i++) {
- cout << expMLintern[i] << " ";
- }
- cout << endl;
-
- cout << "expMLbase[31]:";
- for(int i = 0; i < 31; i++) {
- cout << i << ":" << expMLbase[i] << endl;
- }
-
- cout << "expint11[6][6][16]:";
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 6; j++) {
- for(int k = 0; k < 16; k++) {
- cout << expint11[i][j][k] << " ";
- }
- cout << endl;
- }
- cout << endl;
- }
-
- cout << "expint21[6][6][16][4]:" << endl;
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 6; j++) {
- for(int k = 0; k < 16; k++) {
- for(int l = 0; l < 4; l++) {
- cout << expint21[i][j][k][l] << " ";
- }
- cout << endl;
- }
- cout << endl;
- }
- cout << endl;
- }
-
-
- cout << "expint22[6][6][16][16]:" << endl;
- for(int i = 0; i < 6; i++) {
- for(int j = 0; j < 6; j++) {
- for(int k = 0; k < 16; k++) {
- for(int l = 0; l < 16; l++) {
- cout << expint22[i][j][k][l] << " ";
- }
- cout << endl;
- }
- cout << endl;
- }
- cout << endl;
- }
-
-}
-
-}