// Program 8 - create a pseudoknot using distance bounds and refinement
molecule m;
float 	xyz[ dynamic ],f[ dynamic ],v[ dynamic ];
bounds 	b;
int	i, seqlen, ier;
float	fret;
string	seq, opt;

seq = " aagagagguggagggacuggcccgaugaacccggcaaccagaauggugccauuccugcagcgguuucgcugaaagaugagagauucuuguagucucuucuuuuagcg ";

seqlen = length(seq);

m = link_na("1", seq, "rna.amber94.rlb", "rna", "35");

allocate xyz[ 4*m.natoms ];
allocate f[ 4*m.natoms ];
allocate v[ 4*m.natoms ];

b = newbounds(m, "");

for ( i = 1; i <= seqlen; i = i + 1) {
	useboundsfrom(b, m, sprintf("1:%d:??,H?[^'T]", i), m,
		sprintf("1:%d:??,H?[^'T]", i), 0.0 );
}

setboundsfromdb(b, m, "1:2:", "1:3:", "arna.stack.db",     1.0);
setboundsfromdb(b, m, "1:3:", "1:4:", "arna.stack.db",     1.0);

setboundsfromdb(b, m, "1:7:", "1:8:", "arna.stack.db",     1.0);
setboundsfromdb(b, m, "1:8:", "1:9:", "arna.stack.db",     1.0);
setboundsfromdb(b, m, "1:9:", "1:10:", "arna.stack.db",    1.0);
setboundsfromdb(b, m, "1:13:", "1:14:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:14:", "1:15:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:16:", "1:17:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:18:", "1:19:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:20:", "1:21:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:21:", "1:22:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:28:", "1:29:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:29:", "1:30:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:30:", "1:31:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:31:", "1:32:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:32:", "1:33:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:34:", "1:35:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:35:", "1:36:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:36:", "1:37:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:37:", "1:38:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:40:", "1:41:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:41:", "1:42:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:43:", "1:44:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:44:", "1:45:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:47:", "1:48:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:48:", "1:49:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:49:", "1:50:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:50:", "1:51:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:53:", "1:54:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:58:", "1:59:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:59:", "1:60:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:62:", "1:63:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:63:", "1:64:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:64:", "1:65:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:65:", "1:66:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:67:", "1:68:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:68:", "1:69:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:69:", "1:70:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:71:", "1:72:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:73:", "1:74:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:74:", "1:75:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:75:", "1:76:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:76:", "1:77:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:77:", "1:78:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:78:", "1:79:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:79:", "1:80:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:81:", "1:82:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:82:", "1:83:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:83:", "1:84:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:85:", "1:86:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:86:", "1:87:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:87:", "1:88:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:94:", "1:95:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:95:", "1:96:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:96:", "1:97:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:97:", "1:98:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:98:", "1:99:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:99:", "1:100:", "arna.stack.db",  1.0);
setboundsfromdb(b, m, "1:100:", "1:101:", "arna.stack.db", 1.0);

setboundsfromdb(b, m, "1:103:", "1:104:", "arna.stack.db", 1.0);
setboundsfromdb(b, m, "1:104:", "1:105:", "arna.stack.db", 1.0);
setboundsfromdb(b, m, "1:105:", "1:106:", "arna.stack.db", 1.0);

setboundsfromdb(b, m, "1:5:", "1:34:", "arna.stack.db",    1.0);
setboundsfromdb(b, m, "1:12:", "1:25:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:12:", "1:28:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:13:", "1:24:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:25:", "1:27:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:41:", "1:47:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:54:", "1:57:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:62:", "1:74:", "arna.stack.db",   1.0);
setboundsfromdb(b, m, "1:65:", "1:71:", "arna.stack.db",   1.0);

setboundsfromdb(b, m, "1:5:", "1:33:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:7:", "1:32:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:8:", "1:31:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:9:", "1:30:", "arna.basepair.db", 1.0);

setboundsfromdb(b, m, "1:10:", "1:28:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:10:", "1:29:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:11:", "1:28:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:12:", "1:24:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:12:", "1:27:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:13:", "1:23:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:13:", "1:28:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:14:", "1:22:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:15:", "1:21:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:15:", "1:54:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:16:", "1:77:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:17:", "1:60:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:18:", "1:59:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:19:", "1:58:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:20:", "1:57:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:22:", "1:53:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:34:", "1:52:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:35:", "1:51:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:36:", "1:50:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:39:", "1:48:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:40:", "1:47:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:41:", "1:46:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:42:", "1:45:", "arna.basepair.db", 1.0);

setboundsfromdb(b, m, "1:61:", "1:74:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:62:", "1:73:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:63:", "1:72:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:64:", "1:71:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:65:", "1:70:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:66:", "1:69:", "arna.basepair.db", 1.0);

setboundsfromdb(b, m, "1:76:", "1:106:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:77:", "1:105:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:78:", "1:104:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:79:", "1:103:", "arna.basepair.db", 1.0);

setboundsfromdb(b, m, "1:81:", "1:101:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:82:", "1:100:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:83:", "1:99:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:84:", "1:98:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:85:", "1:97:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:86:", "1:96:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:87:", "1:95:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:88:", "1:94:", "arna.basepair.db", 1.0);
setboundsfromdb(b, m, "1:89:", "1:93:", "arna.basepair.db", 1.0);

tsmooth( b, 0.0005 );

opt = "seed=571, gdist=0, ntpr=50, k4d=2.0, randpair=5";
dg_options( b, opt );
embed(b, xyz );

ier = conjgrad( xyz, 4*m.natoms, fret, db_viol, 0.1, 10., 500 );

for ( i = 3000; i > 2800; i = i - 100 ){
	dg_options( b, "ntpr=1000, k4d=0.2" );
	mm_options( "ntpr_md=50, zerov=1, temp0=" + sprintf("%d.", i) );
	md( 4*m.natoms, 1000, xyz, f, v, db_viol );

	dg_options( b, "ntpr=1000, k4d=4.0" );
	mm_options( "zerov=0, temp0=0., tautp=0.3" );
	md( 4*m.natoms, 8000, xyz, f, v, db_viol );
}

setmol_from_xyzw( m, NULL, xyz );
putpdb( "rna_pseudoknot.pdb", m );