RNAlib-2.2.0-RC3
multibranch_loops.h
Go to the documentation of this file.
1 #ifndef VIENNA_RNA_PACKAGE_MULTIBRANCH_LOOPS_H
2 #define VIENNA_RNA_PACKAGE_MULTIBRANCH_LOOPS_H
3 
4 #include <stdio.h>
5 #include <stdlib.h>
6 #include <math.h>
7 #include <ctype.h>
8 #include <string.h>
9 #include <ViennaRNA/fold_vars.h>
10 #include <ViennaRNA/energy_par.h>
11 #include <ViennaRNA/params.h>
12 #include <ViennaRNA/constraints.h>
14 #include <ViennaRNA/gquad.h>
15 
16 #ifdef __GNUC__
17 # define INLINE inline
18 #else
19 # define INLINE
20 #endif
21 
22 #ifdef ON_SAME_STRAND
23 #undef ON_SAME_STRAND
24 #endif
25 
26 #define ON_SAME_STRAND(I,J,C) (((I)>=(C))||((J)<(C)))
27 
56 INLINE PRIVATE int E_MLstem( int type,
57  int si1,
58  int sj1,
59  vrna_param_t *P);
60 
67 INLINE PRIVATE double exp_E_MLstem(int type,
68  int si1,
69  int sj1,
70  vrna_exp_param_t *P);
71 
72 
73 
79 INLINE PRIVATE int E_mb_loop_stack(int i, int j, vrna_fold_compound *vc);
80 
96 INLINE PRIVATE int
98  int *i,
99  int *j,
100  int *k,
101  int en,
102  int *component1,
103  int *component2);
104 
105 /*
106 #################################
107 # BEGIN OF FUNCTION DEFINITIONS #
108 #################################
109 */
110 
111 
112 INLINE PRIVATE int
113 E_mb_loop_fast( int i,
114  int j,
115  vrna_fold_compound *vc,
116  int *dmli1,
117  int *dmli2){
118 
119  int decomp, en;
120  unsigned char type, tt;
121  short S_i1, S_j1;
122 
123  int cp = vc->cutpoint;
124  char *ptype = vc->ptype;
125  short *S = vc->sequence_encoding;
126  int *indx = vc->jindx;
127  char *hc = vc->hc->matrix;
128  int *hc_up = vc->hc->up_ml;
129  vrna_sc_t *sc = vc->sc;
130  int *fc = vc->matrices->fc;
131  vrna_param_t *P = vc->params;
132 
133  int ij = indx[j] + i;
134  int hc_decompose = hc[ij];
135  int e = INF;
136  int dangle_model = P->model_details.dangles;
137  int *rtype = &(P->model_details.rtype[0]);
138 
139  type = (unsigned char)ptype[ij];
140 
141  if(cp < 0){
142  S_i1 = S[i+1];
143  S_j1 = S[j-1];
144  } else {
145  S_i1 = ON_SAME_STRAND(i, i + 1, cp) ? S[i+1] : -1;
146  S_j1 = ON_SAME_STRAND(j - 1, j, cp) ? S[j-1] : -1;
147  }
148 
149  if(hc_decompose & VRNA_CONSTRAINT_CONTEXT_MB_LOOP){
150  if((S_i1 >= 0) && (S_j1 >= 0)){ /* regular multi branch loop */
151  decomp = dmli1[j-1];
152  tt = rtype[type];
153  switch(dangle_model){
154  /* no dangles */
155  case 0: if(decomp != INF){
156  decomp += E_MLstem(tt, -1, -1, P);
157  if(sc){
158  if(sc->en_basepair)
159  decomp += sc->en_basepair[ij];
160  }
161  }
162  break;
163 
164  /* double dangles */
165  case 2: decomp += E_MLstem(tt, S_j1, S_i1, P);
166  if(sc){
167  if(sc->en_basepair)
168  decomp += sc->en_basepair[ij];
169  }
170  break;
171 
172  /* normal dangles, aka dangles = 1 || 3 */
173  default: if(decomp != INF){
174  decomp += E_MLstem(tt, -1, -1, P);
175  if(sc){
176  if(sc->en_basepair)
177  decomp += sc->en_basepair[ij];
178  }
179  }
180  if(hc_up[i+1]){
181  if(dmli2[j-1] != INF){
182  en = dmli2[j-1] + E_MLstem(tt, -1, S_i1, P) + P->MLbase;
183  if(sc){
184  if(sc->free_energies)
185  en += sc->free_energies[i+1][1];
186 
187  if(sc->en_basepair)
188  en += sc->en_basepair[ij];
189  }
190  decomp = MIN2(decomp, en);
191  }
192  }
193  if(hc_up[j-1] && hc_up[i+1]){
194  if(dmli2[j-2] != INF){
195  en = dmli2[j-2] + E_MLstem(tt, S_j1, S_i1, P) + 2*P->MLbase;
196  if(sc){
197  if(sc->free_energies)
198  en += sc->free_energies[i+1][1]
199  + sc->free_energies[j-1][1];
200 
201  if(sc->en_basepair)
202  en += sc->en_basepair[ij];
203  }
204  decomp = MIN2(decomp, en);
205  }
206  }
207  if(hc_up[j-1]){
208  if(dmli1[j-2] != INF){
209  en = dmli1[j-2] + E_MLstem(tt, S_j1, -1, P) + P->MLbase;
210  if(sc){
211  if(sc->free_energies)
212  en += sc->free_energies[j-1][1];
213 
214  if(sc->en_basepair)
215  en += sc->en_basepair[ij];
216  }
217  decomp = MIN2(decomp, en);
218  }
219  }
220  break;
221  }
222  if(decomp != INF)
223  e = decomp + P->MLclosing;
224  }
225 
226  if(!ON_SAME_STRAND(i, j, cp)){ /* multibrach like cofold structure with cut somewhere between i and j */
227  decomp = fc[i+1] + fc[j-1];
228  tt = rtype[type];
229  switch(dangle_model){
230  case 0: decomp += E_ExtLoop(tt, -1, -1, P);
231  break;
232  case 2: decomp += E_ExtLoop(tt, S_j1, S_i1, P);
233  break;
234  default: decomp += E_ExtLoop(tt, -1, -1, P);
235  if((hc_up[i+1]) && (hc_up[j-1])){
236  en = fc[i+2] + fc[j-2] + E_ExtLoop(tt, S_j1, S_i1, P);
237  decomp = MIN2(decomp, en);
238  }
239  if(hc_up[i+1]){
240  en = fc[i+2] + fc[j-1] + E_ExtLoop(tt, -1, S_i1, P);
241  decomp = MIN2(decomp, en);
242  }
243  if(hc_up[j-1]){
244  en = fc[i+1] + fc[j-2] + E_ExtLoop(tt, S_j1, -1, P);
245  decomp = MIN2(decomp, en);
246  }
247  break;
248  }
249  e = MIN2(e, decomp);
250  }
251  }
252  return e;
253 }
254 
255 INLINE PRIVATE int
257  int j,
258  vrna_fold_compound *vc){
259 
260  int e, decomp, en, i1k, k1j1, ij, k;
261  unsigned char type, type_2;
262 
263  int *indx = vc->jindx;
264  char *hc = vc->hc->matrix;
265  int *c = vc->matrices->c;
266  int *fML = vc->matrices->fML;
267  vrna_param_t *P = vc->params;
268  vrna_md_t *md = &(P->model_details);
269  int turn = md->min_loop_size;
270  char *ptype = vc->ptype;
271  int *rtype = &(md->rtype[0]);
272  vrna_sc_t *sc = vc->sc;
273 
274  e = INF;
275  ij = indx[j] + i;
276  type = ptype[ij];
277 
278  if(hc[ij] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP){
279  decomp = INF;
280  k1j1 = indx[j-1] + i + 2 + turn + 1;
281  for (k = i+2+turn; k < j-2-turn; k++, k1j1++){
282  i1k = indx[k] + i + 1;
284  type_2 = rtype[(unsigned char)ptype[i1k]];
285  en = c[i1k]+P->stack[type][type_2]+fML[k1j1];
286  decomp = MIN2(decomp, en);
287  }
288  if(hc[k1j1] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
289  type_2 = rtype[(unsigned char)ptype[k1j1]];
290  en = c[k1j1]+P->stack[type][type_2]+fML[i1k];
291  decomp = MIN2(decomp, en);
292  }
293  }
294  /* no TermAU penalty if coax stack */
295  decomp += 2*P->MLintern[1] + P->MLclosing;
296  if(sc){
297  if(sc->en_basepair)
298  decomp += sc->en_basepair[ij];
299  }
300  e = decomp;
301 
302  }
303  return e;
304 }
305 
306 INLINE PRIVATE int
307 E_ml_rightmost_stem(int i,
308  int j,
309  vrna_fold_compound *vc){
310 
311  int en;
312  vrna_param_t *P = vc->params;
313  int length = vc->length;
314  short *S = vc->sequence_encoding;
315  int *indx = vc->jindx;
316  char *hc = vc->hc->matrix;
317  int *hc_up = vc->hc->up_ml;
318  vrna_sc_t *sc = vc->sc;
319  int *c = vc->matrices->c;
320  int *fm = (P->model_details.uniq_ML) ? vc->matrices->fM1 : vc->matrices->fML;
321  int *ggg = vc->matrices->ggg;
322  int ij = indx[j] + i;
323  int type = vc->ptype[ij];
324  int hc_decompose = hc[ij];
325  int dangle_model = P->model_details.dangles;
326  int with_gquad = P->model_details.gquad;
327  int cp = vc->cutpoint;
328  int e = INF;
329 
330  if(ON_SAME_STRAND(i - 1, i, cp)){
331  if(ON_SAME_STRAND(j, j + 1, cp)){
332  if(hc_decompose & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
333  e = c[ij];
334  if(e != INF){
335  switch(dangle_model){
336  case 2: e += E_MLstem(type, (i==1) ? S[length] : S[i-1], S[j+1], P);
337  break;
338  default: e += E_MLstem(type, -1, -1, P);
339  break;
340  }
341  }
342  }
343 
344  if(with_gquad)
345  if(ON_SAME_STRAND(i, j, cp)){
346  en = ggg[ij] + E_MLstem(0, -1, -1, P);
347  e = MIN2(e, en);
348  }
349 
350  }
351 
352  if(ON_SAME_STRAND(j - 1, j, cp))
353  if(hc_up[j]){
354  if(fm[indx[j - 1] + i] != INF){
355  en = fm[indx[j - 1] + i] + P->MLbase;
356  if(sc)
357  if(sc->free_energies)
358  en += sc->free_energies[j][1];
359 
360  e = MIN2(e, en);
361  }
362  }
363  }
364  return e;
365 }
366 
367 INLINE PRIVATE int
368 E_ml_stems_fast(int i,
369  int j,
370  vrna_fold_compound *vc,
371  int *fmi,
372  int *dmli){
373 
374  int k, en, decomp, mm5, mm3, type_2, k1j, stop;
375 
376  int length = (int)vc->length;
377  char *ptype = vc->ptype;
378  short *S = vc->sequence_encoding;
379  int *indx = vc->jindx;
380  char *hc = vc->hc->matrix;
381  int *hc_up = vc->hc->up_ml;
382  vrna_sc_t *sc = vc->sc;
383  int *c = vc->matrices->c;
384  int *fm = vc->matrices->fML;
385  vrna_param_t *P = vc->params;
386  int ij = indx[j] + i;
387  int dangle_model = P->model_details.dangles;
388  int turn = P->model_details.min_loop_size;
389  int type = ptype[ij];
390  int *rtype = &(P->model_details.rtype[0]);
391  int circular = P->model_details.circ;
392  int cp = vc->cutpoint;
393  int e = INF;
394 
395  /* extension with one unpaired nucleotide at the right (3' site)
396  or full branch of (i,j)
397  */
398  e = E_ml_rightmost_stem(i,j,vc);
399 
400  /* extension with one unpaired nucleotide at 5' site
401  and all other variants which are needed for odd
402  dangle models
403  */
404  if(ON_SAME_STRAND(i - 1, i, cp)){
405  switch(dangle_model){
406  /* no dangles */
407  case 0: /* fall through */
408 
409  /* double dangles */
410  case 2: if(ON_SAME_STRAND(i, i + 1, cp))
411  if(hc_up[i]){
412  if(fm[ij + 1] != INF){
413  en = fm[ij + 1] + P->MLbase;
414  if(sc)
415  if(sc->free_energies)
416  en += sc->free_energies[i][1];
417  e = MIN2(e, en);
418  }
419  }
420  break;
421 
422  /* normal dangles, aka dangle_model = 1 || 3 */
423  default: mm5 = ((i>1) || circular) ? S[i] : -1;
424  mm3 = ((j<length) || circular) ? S[j] : -1;
425  if(ON_SAME_STRAND(i, i + 1, cp))
426  if(hc_up[i]){
427  if(fm[ij+1] != INF){
428  en = fm[ij+1] + P->MLbase;
429  if(sc)
430  if(sc->free_energies)
431  en += sc->free_energies[i][1];
432  e = MIN2(e, en);
433  }
434  if(hc[ij+1] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
435  if(c[ij+1] != INF){
436  type = ptype[ij+1];
437  en = c[ij+1] + E_MLstem(type, mm5, -1, P) + P->MLbase;
438  if(sc)
439  if(sc->free_energies)
440  en += sc->free_energies[i][1];
441  e = MIN2(e, en);
442  }
443  }
444  }
445 
446  if(ON_SAME_STRAND(j - 1, j, cp))
447  if(hc_up[j]){
448  if(hc[indx[j-1]+i] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
449  if(c[indx[j-1]+i] != INF){
450  type = ptype[indx[j-1]+i];
451  en = c[indx[j-1]+i] + E_MLstem(type, -1, mm3, P) + P->MLbase;
452  if(sc)
453  if(sc->free_energies)
454  en += sc->free_energies[j][1];
455  e = MIN2(e, en);
456  }
457  }
458  }
459 
460  if(ON_SAME_STRAND(j - 1, j, cp) && ON_SAME_STRAND(i, i + 1, cp))
461  if(hc[indx[j-1]+i+1] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
462  if(hc_up[i] && hc_up[j]){
463  if(c[indx[j-1]+i+1] != INF){
464  type = ptype[indx[j-1]+i+1];
465  en = c[indx[j-1]+i+1] + E_MLstem(type, mm5, mm3, P) + 2*P->MLbase;
466  if(sc)
467  if(sc->free_energies)
468  en += sc->free_energies[j][1] + sc->free_energies[i][1];
469  e = MIN2(e, en);
470  }
471  }
472  }
473  break;
474  }
475  }
476 
477  /* modular decomposition -------------------------------*/
478  k1j = indx[j] + i + turn + 2;
479  stop = (cp > 0) ? (cp - 1) : (j - 2 - turn);
480  for (decomp = INF, k = i + 1 + turn; k <= stop; k++, k1j++){
481  if((fmi[k] != INF ) && (fm[k1j] != INF)){
482  en = fmi[k] + fm[k1j];
483  decomp = MIN2(decomp, en);
484  }
485  }
486  k++; k1j++;
487  for (;k <= j - 2 - turn; k++, k1j++){
488  if((fmi[k] != INF) && (fm[k1j] != INF)){
489  en = fmi[k] + fm[k1j];
490  decomp = MIN2(decomp, en);
491  }
492  }
493 
494  dmli[j] = decomp; /* store for use in fast ML decompositon */
495  e = MIN2(e, decomp);
496 
497  /* coaxial stacking */
498  if (dangle_model==3) {
499  /* additional ML decomposition as two coaxially stacked helices */
500  int ik;
501  k1j = indx[j]+i+turn+2;
502  for (decomp = INF, k = i + 1 + turn; k <= stop; k++, k1j++){
503  ik = indx[k]+i;
504  if((hc[ik] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC) && (hc[k1j] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC)){
505  type = rtype[(unsigned char)ptype[ik]];
506  type_2 = rtype[(unsigned char)ptype[k1j]];
507  en = c[ik] + c[k1j] + P->stack[type][type_2];
508  decomp = MIN2(decomp, en);
509  }
510  }
511  k++; k1j++;
512  for (; k <= j-2-turn; k++, k1j++){
513  ik = indx[k]+i;
514  if((hc[ik] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC) && (hc[k1j] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC)){
515  type = rtype[(unsigned char)ptype[ik]];
516  type_2 = rtype[(unsigned char)ptype[k1j]];
517  en = c[ik] + c[k1j] + P->stack[type][type_2];
518  decomp = MIN2(decomp, en);
519  }
520  }
521 
522  decomp += 2*P->MLintern[1]; /* no TermAU penalty if coax stack */
523 #if 0
524  /* This is needed for Y shaped ML loops with coax stacking of
525  interior pairts, but backtracking will fail if activated */
526  DMLi[j] = MIN2(DMLi[j], decomp);
527  DMLi[j] = MIN2(DMLi[j], DMLi[j-1]+P->MLbase);
528  DMLi[j] = MIN2(DMLi[j], DMLi1[j]+P->MLbase);
529  new_fML = MIN2(new_fML, DMLi[j]);
530 #endif
531  e = MIN2(e, decomp);
532  }
533 
534  fmi[j] = e;
535 
536  return e;
537 }
538 
539 
540 
541 INLINE PRIVATE int E_MLstem(int type, int si1, int sj1, vrna_param_t *P){
542  int energy = 0;
543  if(si1 >= 0 && sj1 >= 0){
544  energy += P->mismatchM[type][si1][sj1];
545  }
546  else if (si1 >= 0){
547  energy += P->dangle5[type][si1];
548  }
549  else if (sj1 >= 0){
550  energy += P->dangle3[type][sj1];
551  }
552 
553  if(type > 2)
554  energy += P->TerminalAU;
555 
556  energy += P->MLintern[type];
557 
558  return energy;
559 }
560 
561 
562 
563 INLINE PRIVATE double
564 exp_E_MLstem( int type,
565  int si1,
566  int sj1,
567  vrna_exp_param_t *P){
568 
569  double energy = 1.0;
570  if(si1 >= 0 && sj1 >= 0){
571  energy = P->expmismatchM[type][si1][sj1];
572  }
573  else if(si1 >= 0){
574  energy = P->expdangle5[type][si1];
575  }
576  else if(sj1 >= 0){
577  energy = P->expdangle3[type][sj1];
578  }
579 
580  if(type > 2)
581  energy *= P->expTermAU;
582 
583  energy *= P->expMLintern[type];
584  return energy;
585 }
586 
587 /*
588 #################################
589 # Backtracking functions below #
590 #################################
591 */
592 
593 INLINE PRIVATE int
594 vrna_BT_mb_loop_fake( vrna_fold_compound *vc,
595  int *u,
596  int *i,
597  int *j,
598  bondT *bp_stack,
599  int *stack_count){
600 
601  int length, ii, jj, k, en, cp, fij, fi, *my_c, *my_fc, *my_ggg, *idx, with_gquad, dangle_model, turn;
602  unsigned char type;
603  char *ptype;
604  short mm5, mm3, *S1;
605  vrna_param_t *P;
606  vrna_md_t *md;
607  vrna_hc_t *hc;
608  vrna_sc_t *sc;
609 
610  cp = vc->cutpoint;
611  length = vc->length;
612  P = vc->params;
613  md = &(P->model_details);
614  hc = vc->hc;
615  sc = vc->sc;
616  S1 = vc->sequence_encoding;
617  ptype = vc->ptype;
618  idx = vc->jindx;
619  my_c = vc->matrices->c;
620  my_fc = vc->matrices->fc;
621  my_ggg = vc->matrices->ggg;
622  turn = md->min_loop_size;
623  with_gquad = md->gquad;
624  dangle_model = md->dangles;
625 
626  ii = *i;
627  jj = *j;
628 
629  if(ii < cp){ /* 'lower' part (fc[i<cut,j=cut-1]) */
630 
631  /* nibble off unpaired 5' bases */
632  do{
633  fij = my_fc[ii];
634  fi = (hc->up_ext[ii]) ? my_fc[ii + 1] : INF;
635 
636  if(sc){
637  if(sc->free_energies)
638  fi += sc->free_energies[ii][1];
639  }
640 
641  if(++ii == jj)
642  break;
643 
644  } while (fij == fi);
645  ii--;
646 
647  if (jj < ii + turn + 2){ /* no more pairs */
648  *u = *i = *j = -1;
649  return 1;
650  }
651 
652  mm5 = (ii > 1 && ON_SAME_STRAND(ii - 1, ii, cp)) ? S1[ii - 1] : -1;
653 
654  /* i or i+1 is paired. Find pairing partner */
655  switch(dangle_model){
656  case 0: for(k = ii + turn + 1; k <= jj; k++){
657  if(hc->matrix[idx[k] + ii] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
658  type = (unsigned char)ptype[idx[k] + ii];
659  if(fij == my_fc[k + 1] + my_c[idx[k] + ii] + E_ExtLoop(type, -1, -1, P)){
660  bp_stack[++(*stack_count)].i = ii;
661  bp_stack[(*stack_count)].j = k;
662  *u = k + 1;
663  *i = ii;
664  *j = k;
665  return 1;
666  }
667  }
668 
669  if (with_gquad){
670  if(fij == my_fc[k + 1] + my_ggg[idx[k] + ii]){
671  *u = k + 1;
672  *i = *j = -1;
673  vrna_BT_gquad_mfe(vc, ii, k, bp_stack, stack_count);
674  return 1;
675  }
676  }
677  }
678  break;
679 
680  case 2: for(k = ii + turn + 1; k <= jj; k++){
681  if(hc->matrix[idx[k] + ii] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
682  mm3 = ON_SAME_STRAND(k, k + 1, cp) ? S1[k + 1] : -1;
683  type = (unsigned char)ptype[idx[k] + ii];
684  if(fij == my_fc[k + 1] + my_c[idx[k] + ii] + E_ExtLoop(type, mm5, mm3, P)){
685  bp_stack[++(*stack_count)].i = ii;
686  bp_stack[(*stack_count)].j = k;
687  *u = k + 1;
688  *i = ii;
689  *j = k;
690  return 1;
691  }
692  }
693 
694  if (with_gquad){
695  if(fij == my_fc[k + 1] + my_ggg[idx[k] + ii]){
696  *u = k + 1;
697  *i = *j = -1;
698  vrna_BT_gquad_mfe(vc, ii, k, bp_stack, stack_count);
699  return 1;
700  }
701  }
702  }
703  break;
704 
705  default: for(k = ii + turn + 1; k <= jj; k++){
706  if(hc->matrix[idx[k] + ii] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
707  type = (unsigned char)ptype[idx[k] + ii];
708  if(fij == my_fc[k + 1] + my_c[idx[k] + ii] + E_ExtLoop(type, -1, -1, P)){
709  bp_stack[++(*stack_count)].i = ii;
710  bp_stack[(*stack_count)].j = k;
711  *u = k + 1;
712  *i = ii;
713  *j = k;
714  return 1;
715  }
716  if(hc->up_ext[k + 1]){
717  mm3 = ON_SAME_STRAND(k, k + 1, cp) ? S1[k + 1] : -1;
718  en = my_c[idx[k] + ii];
719  if(sc)
720  if(sc->free_energies)
721  en += sc->free_energies[k + 1][1];
722 
723  if(fij == my_fc[k + 2] + en + E_ExtLoop(type, -1, mm3, P)){
724  bp_stack[++(*stack_count)].i = ii;
725  bp_stack[(*stack_count)].j = k;
726  *u = k + 2;
727  *i = ii;
728  *j = k;
729  return 1;
730  }
731  }
732  }
733 
734  if (with_gquad){
735  if(fij == my_fc[k + 1] + my_ggg[idx[k] + ii]){
736  *u = k + 1;
737  *i = *j = -1;
738  vrna_BT_gquad_mfe(vc, ii, k, bp_stack, stack_count);
739  return 1;
740  }
741  }
742 
743  if(hc->matrix[idx[k] + ii + 1] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
744  if(hc->up_ext[ii]){
745  mm5 = ON_SAME_STRAND(ii, ii + 1,cp) ? S1[ii] : -1;
746  mm3 = ON_SAME_STRAND(k, k + 1,cp) ? S1[k + 1] : -1;
747  type = ptype[idx[k] + ii + 1];
748  en = my_c[idx[k] + ii + 1];
749  if(sc)
750  if(sc->free_energies)
751  en += sc->free_energies[ii][1];
752 
753  if(fij == en + my_fc[k + 1] + E_ExtLoop(type, mm5, -1, P)){
754  bp_stack[++(*stack_count)].i = ii + 1;
755  bp_stack[(*stack_count)].j = k;
756  *u = k + 1;
757  *i = ii + 1;
758  *j = k;
759  return 1;
760  }
761 
762  if(k < jj){
763  if(hc->up_ext[k + 1]){
764  if(sc)
765  if(sc->free_energies)
766  en += sc->free_energies[k + 1][1];
767 
768  if(fij == en + my_fc[k + 2] + E_ExtLoop(type, mm5, mm3, P)){
769  bp_stack[++(*stack_count)].i = ii + 1;
770  bp_stack[(*stack_count)].j = k;
771  *u = k + 2;
772  *i = ii + 1;
773  *j = k;
774  return 1;
775  }
776  }
777  }
778  }
779  }
780  }
781  break;
782  }
783  } else { /* 'upper' part (fc[i=cut,j>cut]) */
784 
785  /* nibble off unpaired 3' bases */
786  do{
787  fij = my_fc[jj];
788  fi = (hc->up_ext[jj]) ? my_fc[jj - 1] : INF;
789 
790  if(sc){
791  if(sc->free_energies)
792  fi += sc->free_energies[jj][1];
793  }
794 
795  if(--jj == ii)
796  break;
797 
798  } while (fij == fi);
799  jj++;
800 
801  if (jj < ii + turn + 2){ /* no more pairs */
802  *u = *i = *j = -1;
803  return 1;
804  }
805 
806  /* j or j-1 is paired. Find pairing partner */
807  mm3 = ((jj < length) && ON_SAME_STRAND(jj, jj + 1, cp)) ? S1[jj + 1] : -1;
808  switch(dangle_model){
809  case 0: for (k = jj - turn - 1; k >= ii; k--){
810  if(with_gquad){
811  if(fij == my_fc[k - 1] + my_ggg[idx[jj] + k]){
812  *u = k - 1;
813  *i = *j = -1;
814  vrna_BT_gquad_mfe(vc, k, jj, bp_stack, stack_count);
815  return 1;
816  }
817  }
818 
819  if(hc->matrix[idx[jj] + k] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
820  type = (unsigned char)ptype[idx[jj] + k];
821  en = my_c[idx[jj] + k];
822  if(!ON_SAME_STRAND(k, jj, cp))
823  en += P->DuplexInit;
824 
825  if(fij == my_fc[k - 1] + en + E_ExtLoop(type, -1, -1, P)){
826  bp_stack[++(*stack_count)].i = k;
827  bp_stack[(*stack_count)].j = jj;
828  *u = k - 1;
829  *i = k;
830  *j = jj;
831  return 1;
832  }
833  }
834  }
835  break;
836  case 2: for(k = jj - turn - 1; k >= ii; k--){
837  if(with_gquad){
838  if(fij == my_fc[k - 1] + my_ggg[idx[jj] + k]){
839  *u = k - 1;
840  *i = *j = -1;
841  vrna_BT_gquad_mfe(vc, k, jj, bp_stack, stack_count);
842  return 1;
843  }
844  }
845 
846  if(hc->matrix[idx[jj] + k] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
847  mm5 = ((k > 1) && ON_SAME_STRAND(k - 1, k, cp)) ? S1[k - 1] : -1;
848  type = (unsigned char)ptype[idx[jj] + k];
849  en = my_c[idx[jj] + k];
850  if(!ON_SAME_STRAND(k, jj, cp))
851  en += P->DuplexInit;
852 
853  if(fij == my_fc[k - 1] + en + E_ExtLoop(type, mm5, mm3, P)){
854  bp_stack[++(*stack_count)].i = k;
855  bp_stack[(*stack_count)].j = jj;
856  *u = k - 1;
857  *i = k;
858  *j = jj;
859  return 1;
860  }
861  }
862  }
863  break;
864 
865  default: for(k = jj - turn - 1; k >= ii; k--){
866 
867  if(with_gquad){
868  if(fij == my_fc[k - 1] + my_ggg[idx[jj] + k]){
869  *u = k - 1;
870  *i = *j = -1;
871  vrna_BT_gquad_mfe(vc, k, jj, bp_stack, stack_count);
872  return 1;
873  }
874  }
875 
876  if(hc->matrix[idx[jj] + k] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
877  type = (unsigned char)ptype[idx[jj] + k];
878  en = my_c[idx[jj] + k];
879  if(!ON_SAME_STRAND(k, jj, cp))
880  en += P->DuplexInit;
881 
882  if(fij == my_fc[k - 1] + en + E_ExtLoop(type, -1, -1, P)){
883  bp_stack[++(*stack_count)].i = k;
884  bp_stack[(*stack_count)].j = jj;
885  *u = k - 1;
886  *i = k;
887  *j = jj;
888  return 1;
889  }
890  if(hc->up_ext[k - 1]){
891  if((k > 1) && ON_SAME_STRAND(k - 1, k, cp)){
892  mm5 = S1[k - 1];
893  if(sc)
894  if(sc->free_energies)
895  en += sc->free_energies[k - 1][1];
896 
897  if(fij == my_fc[k - 2] + en + E_ExtLoop(type, mm5, -1, P)){
898  bp_stack[++(*stack_count)].i = k;
899  bp_stack[(*stack_count)].j = jj;
900  *u = k - 2;
901  *i = k;
902  *j = jj;
903  return 1;
904  }
905  }
906  }
907  }
908 
909  if(hc->matrix[idx[jj - 1] + k] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
910  type = (unsigned char)ptype[idx[jj - 1] + k];
911  if(hc->up_ext[jj]){
912  if(ON_SAME_STRAND(jj - 1, jj, cp)){
913  mm3 = S1[jj];
914  en = my_c[idx[jj - 1] + k];
915  if (!ON_SAME_STRAND(k, jj - 1, cp))
916  en += P->DuplexInit; /*???*/
917  if(sc)
918  if(sc->free_energies)
919  en += sc->free_energies[jj][1];
920 
921  if (fij == en + my_fc[k - 1] + E_ExtLoop(type, -1, mm3, P)){
922  bp_stack[++(*stack_count)].i = k;
923  bp_stack[(*stack_count)].j = jj - 1;
924  *u = k - 1;
925  *i = k;
926  *j = jj - 1;
927  return 1;
928  }
929 
930  if(k > ii){
931  if(hc->up_ext[k - 1]){
932  mm5 = ON_SAME_STRAND(k - 1, k, cp) ? S1[k - 1] : -1;
933  if(sc)
934  if(sc->free_energies)
935  en += sc->free_energies[k - 1][1];
936 
937  if (fij == my_fc[k - 2] + en + E_ExtLoop(type, mm5, mm3, P)){
938  bp_stack[++(*stack_count)].i = k;
939  bp_stack[(*stack_count)].j = jj - 1;
940  *u = k - 2;
941  *i = k;
942  *j = jj - 1;
943  return 1;
944  }
945  }
946  }
947  }
948  }
949  }
950  }
951  break;
952  }
953  }
954 
955  return 0;
956 }
957 
958 INLINE PRIVATE int
959 vrna_BT_mb_loop_split(vrna_fold_compound *vc,
960  int *i,
961  int *j,
962  int *k,
963  int *l,
964  int *component1,
965  int *component2,
966  bondT *bp_stack,
967  int *stack_count){
968 
969  int length, cp, ij, ii, jj, fij, fi, u, en, *my_c, *my_fML, *my_ggg, turn, *idx, with_gquad, dangle_model, *rtype;
970  unsigned char type, type_2;
971  char *ptype;
972  short *S1;
973  vrna_param_t *P;
974  vrna_md_t *md;
975  vrna_hc_t *hc;
976  vrna_sc_t *sc;
977 
978  length = vc->length;
979  cp = vc->cutpoint;
980  P = vc->params;
981  md = &(P->model_details);
982  hc = vc->hc;
983  sc = vc->sc;
984  idx = vc->jindx;
985  ptype = vc->ptype;
986  rtype = &(md->rtype[0]);
987  S1 = vc->sequence_encoding;
988 
989  my_c = vc->matrices->c;
990  my_fML = vc->matrices->fML;
991  my_ggg = vc->matrices->ggg;
992  turn = md->min_loop_size;
993  with_gquad = md->gquad;
994  dangle_model = md->dangles;
995 
996  ii = *i;
997  jj = *j;
998 
999  /* nibble off unpaired 3' bases */
1000  do{
1001  fij = my_fML[idx[jj] + ii];
1002  fi = (hc->up_ml[jj]) ? my_fML[idx[jj - 1] + ii] + P->MLbase : INF;
1003 
1004  if(sc){
1005  if(sc->free_energies)
1006  fi += sc->free_energies[jj][1];
1007  }
1008 
1009  if(--jj == 0)
1010  break;
1011 
1012  } while (fij == fi);
1013  jj++;
1014 
1015  /* nibble off unpaired 5' bases */
1016  do{
1017  fij = my_fML[idx[jj] + ii];
1018  fi = (hc->up_ml[ii]) ? my_fML[idx[jj] + ii + 1] + P->MLbase : INF;
1019 
1020  if(sc){
1021  if(sc->free_energies)
1022  fi += sc->free_energies[ii][1];
1023  }
1024 
1025  if(++ii == jj)
1026  break;
1027 
1028  } while (fij == fi);
1029  ii--;
1030 
1031  if (jj < ii + turn + 2){ /* no more pairs */
1032  return 0;
1033  }
1034 
1035  ij = idx[jj] + ii;
1036 
1037  *component1 = *component2 = 1; /* split into two multi loop parts by default */
1038 
1039  /* 1. test for single component */
1040 
1041  if(with_gquad){
1042  if(fij == my_ggg[ij] + E_MLstem(0, -1, -1, P)){
1043  *i = *j = -1;
1044  *k = *l = -1;
1045  vrna_BT_gquad_mfe(vc, ii, jj, bp_stack, stack_count);
1046  return 1;
1047  }
1048  }
1049 
1050  type = (unsigned char)ptype[ij];
1051  en = my_c[ij];
1052 
1053  switch(dangle_model){
1054  case 0: if(hc->matrix[ij] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1055  if(fij == en + E_MLstem(type, -1, -1, P)){
1056  *i = *j = -1;
1057  *k = ii;
1058  *l = jj;
1059  *component2 = 2; /* 2nd part is structure enclosed by base pair */
1060  return 1;
1061  }
1062  }
1063  break;
1064 
1065  case 2: if(hc->matrix[ij] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1066  if(fij == en + E_MLstem(type, S1[ii - 1], S1[jj + 1], P)){
1067  *i = *j = -1;
1068  *k = ii;
1069  *l = jj;
1070  *component2 = 2;
1071  return 1;
1072  }
1073  }
1074  break;
1075 
1076  default: if(hc->matrix[ij] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1077  if(fij == en + E_MLstem(type, -1, -1, P)){
1078  *i = *j = -1;
1079  *k = ii;
1080  *l = jj;
1081  *component2 = 2;
1082  return 1;
1083  }
1084  }
1085  if(hc->matrix[ij + 1] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1086  if(hc->up_ml[ii]){
1087  int tmp_en = fij;
1088  if(sc)
1089  if(sc->free_energies)
1090  tmp_en -= sc->free_energies[ii][1];
1091 
1092  type = (unsigned char)ptype[ij + 1];
1093  if(tmp_en == my_c[ij+1] + E_MLstem(type, S1[ii], -1, P) + P->MLbase){
1094  *i = *j = -1;
1095  *k = ii + 1;
1096  *l = jj;
1097  *component2 = 2;
1098  return 1;
1099  }
1100  }
1101  }
1102  if(hc->matrix[idx[jj - 1] + ii] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1103  if(hc->up_ml[jj]){
1104  int tmp_en = fij;
1105  if(sc)
1106  if(sc->free_energies)
1107  tmp_en -= sc->free_energies[jj][1];
1108 
1109  type = (unsigned char)ptype[idx[jj - 1] + ii];
1110  if(tmp_en == my_c[idx[jj - 1] + ii] + E_MLstem(type, -1, S1[jj], P) + P->MLbase){
1111  *i = *j = -1;
1112  *k = ii;
1113  *l = jj - 1;
1114  *component2 = 2;
1115  return 1;
1116  }
1117  }
1118  }
1119  if(hc->matrix[idx[jj - 1] + ii + 1] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1120  if(hc->up_ml[ii] && hc->up_ml[jj]){
1121  int tmp_en = fij;
1122  if(sc)
1123  if(sc->free_energies)
1124  tmp_en -= sc->free_energies[ii][1] + sc->free_energies[jj][1];
1125 
1126  type = (unsigned char)ptype[idx[jj - 1] + ii + 1];
1127  if(tmp_en == my_c[idx[jj - 1] + ii + 1] + E_MLstem(type, S1[ii], S1[jj], P) + 2 * P->MLbase){
1128  *i = *j = -1;
1129  *k = ii + 1;
1130  *l = jj - 1;
1131  *component2 = 2;
1132  return 1;
1133  }
1134  }
1135  }
1136  break;
1137  }
1138 
1139  /* 2. Test for possible split point */
1140  for(u = ii + 1 + turn; u <= jj - 2 - turn; u++){
1141  if(fij == my_fML[idx[u] + ii] + my_fML[idx[jj] + u + 1]){
1142  *i = ii;
1143  *j = u;
1144  *k = u + 1;
1145  *l = jj;
1146  return 1;
1147  }
1148  }
1149 
1150  /* 3. last chance! Maybe coax stack */
1151  if(dangle_model==3){
1152  int ik, k1j, tmp_en;
1153  for (k1j = idx[jj] + ii + turn + 2, u = ii + 1 + turn; u <= jj - 2 - turn; u++, k1j++) {
1154  ik = idx[u] + ii;
1155  if((hc->matrix[ik] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC) && (hc->matrix[k1j] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC)){
1156  type = rtype[(unsigned char)ptype[ik]];
1157  type_2 = rtype[(unsigned char)ptype[k1j]];
1158  tmp_en = my_c[ik] + my_c[k1j] + P->stack[type][type_2] + 2*P->MLintern[1];
1159  if (fij == tmp_en){
1160  *i = ii;
1161  *j = u;
1162  *k = u + 1;
1163  *l = jj;
1164  *component1 = *component2 = 2;
1165  return 1;
1166  }
1167  }
1168  }
1169  }
1170 
1171  return 0;
1172 }
1173 
1174 
1175 INLINE PRIVATE int
1177  int *i,
1178  int *j,
1179  int *k,
1180  int en,
1181  int *component1,
1182  int *component2){
1183 
1184  int ij, p, q, r, e, cp, *idx, turn, dangle_model, *my_c, *my_fML, *my_fc, *rtype;
1185  unsigned char type, type_2, tt;
1186  char *ptype;
1187  short s5, s3, *S1;
1188  vrna_param_t *P;
1189  vrna_md_t *md;
1190  vrna_hc_t *hc;
1191  vrna_sc_t *sc;
1192 
1193  cp = vc->cutpoint;
1194  idx = vc->jindx;
1195  ij = idx[*j] + *i;
1196  S1 = vc->sequence_encoding;
1197  P = vc->params;
1198  md = &(P->model_details);
1199  hc = vc->hc;
1200  sc = vc->sc;
1201  my_c = vc->matrices->c;
1202  my_fML = vc->matrices->fML;
1203  my_fc = vc->matrices->fc;
1204  turn = md->min_loop_size;
1205  ptype = vc->ptype;
1206  rtype = &(md->rtype[0]);
1207  type = (unsigned char)ptype[ij];
1208  tt = type;
1209  type = rtype[type];
1210  dangle_model = md->dangles;
1211 
1212  p = *i + 1;
1213  q = *j - 1;
1214 
1215  r = q - turn - 1;
1216 
1217  if(hc->matrix[ij] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP){
1218 
1219  /* is it a fake multi-loop? */
1220  if(!ON_SAME_STRAND(*i, *j, cp)){
1221  int ii, jj;
1222  ii = jj = 0;
1223  e = my_fc[p] + my_fc[q];
1224  if(sc)
1225  if(sc->en_basepair)
1226  e += sc->en_basepair[ij];
1227 
1228  s5 = ON_SAME_STRAND(q, *j, cp) ? S1[q] : -1;
1229  s3 = ON_SAME_STRAND(*i, p, cp) ? S1[p] : -1;
1230 
1231  switch(dangle_model){
1232  case 0: if(en == e + E_ExtLoop(type, -1, -1, P)){
1233  ii=p, jj=q;
1234  }
1235  break;
1236  case 2: if(en == e + E_ExtLoop(type, s5, s3, P)){
1237  ii=p, jj=q;
1238  }
1239  break;
1240  default: {
1241  if(en == e + E_ExtLoop(type, -1, -1, P)){
1242  ii=p, jj=q;
1243  break;
1244  }
1245  if(hc->up_ext[p]){
1246  e = my_fc[p + 1] + my_fc[q];
1247  if(sc){
1248  if(sc->free_energies)
1249  e += sc->free_energies[p][1];
1250  if(sc->en_basepair)
1251  e += sc->en_basepair[ij];
1252  }
1253  if(en == e + E_ExtLoop(type, -1, s3, P)){
1254  ii = p + 1; jj = q;
1255  break;
1256  }
1257  }
1258  if(hc->up_ext[q]){
1259  e = my_fc[p] + my_fc[q - 1];
1260  if(sc){
1261  if(sc->free_energies)
1262  e += sc->free_energies[q][1];
1263  if(sc->en_basepair)
1264  e += sc->en_basepair[ij];
1265  }
1266  if(en == e + E_ExtLoop(type, s5, -1, P)){
1267  ii = p; jj = q - 1;
1268  break;
1269  }
1270  }
1271  if((hc->up_ext[q]) && (hc->up_ext[p])){
1272  e = my_fc[p + 1] + my_fc[q - 1];
1273  if(sc){
1274  if(sc->free_energies)
1275  e += sc->free_energies[p][1] + sc->free_energies[q][1];
1276  if(sc->en_basepair)
1277  e += sc->en_basepair[ij];
1278  }
1279  if(en == e + E_ExtLoop(type, s5, s3, P)){
1280  ii = p + 1; jj = q - 1;
1281  break;
1282  }
1283  }
1284  }
1285  break;
1286  }
1287 
1288  if(ii){ /* found a decomposition */
1289  *component1 = 3;
1290  *i = ii;
1291  *k = cp - 1;
1292  *j = jj;
1293  *component2 = 4;
1294  return 1;
1295  }
1296  }
1297 
1298  /* true multi loop? */
1299  *component1 = *component2 = 1; /* both components are MB loop parts by default */
1300 
1301  s5 = ON_SAME_STRAND(q, *j, cp) ? S1[q] : -1;
1302  s3 = ON_SAME_STRAND(*i, p, cp) ? S1[p] : -1;
1303 
1304  switch(dangle_model){
1305  case 0: e = en - E_MLstem(type, -1, -1, P) - P->MLclosing;
1306  if(sc){
1307  if(sc->en_basepair)
1308  e -= sc->en_basepair[ij];
1309  }
1310  for(r = *i + 2 + turn; r < *j - 2 - turn; ++r){
1311  if(e == my_fML[idx[r] + p] + my_fML[idx[q] + r + 1])
1312  break;
1313  }
1314  break;
1315 
1316  case 2: e = en - E_MLstem(type, s5, s3, P) - P->MLclosing;
1317  if(sc){
1318  if(sc->en_basepair)
1319  e -= sc->en_basepair[ij];
1320  }
1321  for(r = p + turn + 1; r < q - turn - 1; ++r){
1322  if(e == my_fML[idx[r] + p] + my_fML[idx[q] + r + 1])
1323  break;
1324  }
1325  break;
1326 
1327  default: e = en - P->MLclosing;
1328  if(sc){
1329  if(sc->en_basepair)
1330  e -= sc->en_basepair[ij];
1331  }
1332  for(r = p + turn + 1; r < q - turn - 1; ++r){
1333  if(e == my_fML[idx[r] + p] + my_fML[idx[q] + r + 1] + E_MLstem(type, -1, -1, P)){
1334  break;
1335  }
1336  if(hc->up_ml[p]){
1337  int tmp_en = e;
1338  if(sc)
1339  if(sc->free_energies)
1340  tmp_en -= sc->free_energies[p][1];
1341 
1342  if(tmp_en == my_fML[idx[r] + p + 1] + my_fML[idx[q] + r + 1] + E_MLstem(type, -1, s3, P) + P->MLbase){
1343  p += 1;
1344  break;
1345  }
1346  }
1347  if(hc->up_ml[q]){
1348  int tmp_en = e;
1349  if(sc)
1350  if(sc->free_energies)
1351  tmp_en -= sc->free_energies[q][1];
1352 
1353  if(tmp_en == my_fML[idx[r] + p] + my_fML[idx[q - 1] + r + 1] + E_MLstem(type, s5, -1, P) + P->MLbase){
1354  q -= 1;
1355  break;
1356  }
1357  }
1358  if((hc->up_ml[p]) && (hc->up_ml[q])){
1359  int tmp_en = e;
1360  if(sc)
1361  if(sc->free_energies)
1362  tmp_en -= sc->free_energies[p][1] + sc->free_energies[q][1];
1363 
1364  if(tmp_en == my_fML[idx[r] + p + 1] + my_fML[idx[q - 1] + r + 1] + E_MLstem(type, s5, s3, P) + 2 * P->MLbase){
1365  p += 1;
1366  q -= 1;
1367  break;
1368  }
1369  }
1370  /* coaxial stacking of (i.j) with (i+1.r) or (r.j-1) */
1371  /* use MLintern[1] since coax stacked pairs don't get TerminalAU */
1372  if(dangle_model == 3){
1373  int tmp_en = e;
1374  if(hc->matrix[idx[r] + p] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1375  type_2 = rtype[(unsigned char)ptype[idx[r] + p]];
1376  tmp_en = my_c[idx[r] + p] + P->stack[tt][type_2] + my_fML[idx[q] + r + 1];
1377  if(e == tmp_en + 2 * P->MLintern[1]){
1378  *component1 = 2;
1379  break;
1380  }
1381  }
1382 
1383  if(hc->matrix[idx[q] + r + 1] & VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC){
1384  type_2 = rtype[(unsigned char)ptype[idx[q] + r + 1]];
1385  tmp_en = my_c[idx[q] + r + 1] + P->stack[tt][type_2] + my_fML[idx[r] + p];
1386  if (e == tmp_en + 2 * P->MLintern[1]){
1387  *component2 = 2;
1388  break;
1389  }
1390  }
1391  }
1392  }
1393  break;
1394  }
1395  }
1396 
1397  if(r <= *j - turn - 3){
1398  *i = p;
1399  *k = r;
1400  *j = q;
1401  return 1;
1402  } else {
1403 #if 0
1404  /* Y shaped ML loops fon't work yet */
1405  if (dangle_model==3) {
1406  d5 = P->dangle5[tt][S1[j-1]];
1407  d3 = P->dangle3[tt][S1[i+1]];
1408  /* (i,j) must close a Y shaped ML loop with coax stacking */
1409  if (cij == fML[indx[j-2]+i+2] + mm + d3 + d5 + P->MLbase + P->MLbase) {
1410  i1 = i+2;
1411  j1 = j-2;
1412  } else if (cij == fML[indx[j-2]+i+1] + mm + d5 + P->MLbase)
1413  j1 = j-2;
1414  else if (cij == fML[indx[j-1]+i+2] + mm + d3 + P->MLbase)
1415  i1 = i+2;
1416  else /* last chance */
1417  if (cij != fML[indx[j-1]+i+1] + mm + P->MLbase)
1418  fprintf(stderr, "backtracking failed in repeat");
1419  /* if we arrive here we can express cij via fML[i1,j1]+dangles */
1420  bt_stack[++s].i = i1;
1421  bt_stack[s].j = j1;
1422  }
1423 #endif
1424  }
1425 
1426  return 0;
1427 }
1428 
1429 
1434 #endif
int * ggg
Energies of g-quadruplexes.
Definition: data_structures.h:432
int uniq_ML
Flag to ensure unique multibranch loop decomposition during folding.
Definition: model.h:186
int * c
Energy array, given that i-j pair.
Definition: data_structures.h:425
struct vrna_hc_t * hc
The hard constraints data structure used for structure prediction.
Definition: data_structures.h:716
int * fM1
Second ML array, only for unique multibrnach loop decomposition.
Definition: data_structures.h:430
int cutpoint
The position of the (cofold) cutpoint within the provided sequence. If there is no cutpoint...
Definition: data_structures.h:712
The hard constraints data structure.
Definition: constraints.h:390
char * matrix
Upper triangular matrix encoding where a base pair or unpaired nucleotide is allowed.
Definition: constraints.h:391
struct vrna_param_t * params
The precomputed free energy contributions for each type of loop.
Definition: data_structures.h:721
int * jindx
DP matrix accessor.
Definition: data_structures.h:725
int * fML
Multi-loop auxiliary energy array.
Definition: data_structures.h:429
int * up_ml
A linear array that holds the number of allowed unpaired nucleotides in a multi branched loop...
Definition: constraints.h:403
int gquad
Include G-quadruplexes in structure prediction.
Definition: model.h:184
int * up_ext
A linear array that holds the number of allowed unpaired nucleotides in an exterior loop...
Definition: constraints.h:394
int ** free_energies
Energy contribution for unpaired sequence stretches.
Definition: constraints.h:414
#define MIN2(A, B)
Get the minimum of two comparable values.
Definition: utils.h:118
vrna_md_t model_details
Model details to be used in the recursions.
Definition: params.h:79
The most basic data structure required by many functions throughout the RNAlib.
Definition: data_structures.h:698
#define VRNA_CONSTRAINT_CONTEXT_MB_LOOP
Hard constraints flag, base pair is enclosed in an interior loop.
Definition: constraints.h:233
PRIVATE int E_ExtLoop(int type, int si1, int sj1, vrna_param_t *P)
Definition: exterior_loops.h:460
int * fc
Energy from i to cutpoint (and vice versa if i>cut)
Definition: data_structures.h:428
#define VRNA_CONSTRAINT_CONTEXT_EXT_LOOP
Hard constraints flag, base pair in the exterior loop.
Definition: constraints.h:201
#define INF
Definition: energy_const.h:16
Energy evaluation of exterior loops for MFE and partition function calculations.
#define VRNA_CONSTRAINT_CONTEXT_MB_LOOP_ENC
Hard constraints flag, base pair is enclosed in a multi branch loop.
Definition: constraints.h:241
vrna_mx_mfe_t * matrices
The MFE DP matrices.
Definition: data_structures.h:718
int min_loop_size
Minimum size of hairpin loops.
Definition: model.h:194
PRIVATE int vrna_BT_mb_loop(vrna_fold_compound *vc, int *i, int *j, int *k, int en, int *component1, int *component2)
Backtrack the decomposition of a multi branch loop closed by .
Definition: multibranch_loops.h:1176
int dangles
Specifies the dangle model used in any energy evaluation (0,1,2 or 3)
Definition: model.h:172
unsigned int length
The length of the sequence (or sequence alignment)
Definition: data_structures.h:711
Various functions related to G-quadruplex computations.
char * ptype
Pair type array.
Definition: data_structures.h:749
Base pair.
Definition: data_structures.h:73
The datastructure that contains temperature scaled energy parameters.
Definition: params.h:41
The soft constraints data structure.
Definition: constraints.h:413
int circ
Assume RNA to be circular instead of linear.
Definition: model.h:183
The data structure that contains the complete model details used throughout the calculations.
Definition: model.h:169
PRIVATE int E_mb_loop_stack(int i, int j, vrna_fold_compound *vc)
Evaluate energy of a multi branch helices stacking onto closing pair (i,j)
Definition: multibranch_loops.h:256
Here all all declarations of the global variables used throughout RNAlib.
short * sequence_encoding
Numerical encoding of the input sequence.
Definition: data_structures.h:744
struct vrna_sc_t * sc
The soft constraints for usage in structure prediction and evaluation.
Definition: data_structures.h:763
int * en_basepair
Energy contribution for base pairs.
Definition: constraints.h:415
The datastructure that contains temperature scaled Boltzmann weights of the energy parameters...
Definition: params.h:86