Krnl_IntraE.cl 12.2 KB
Newer Older
1
2
3
4
5
6
// sqrt7 ////https://www.codeproject.com/Articles/69941/Best-Square-Root-Method-Algorithm-Function-Precisi
float sqrt_custom(const float x) 
{ 	//uint i = as_uint(x);	
	uint i = *(uint*) &x;    	
	i  += 127 << 23;	// adjust bias   	
	i >>= 1; 		// approximation of square root 	
7
	return as_float(i);	// return *(float*) &i; 
8
9
}  

10
// --------------------------------------------------------------------------
11
12
// IntraE calculates the intramolecular energy of a set of atomic ligand 
// contributor-pairs.
13
14
// Originally from: processligand.c
// --------------------------------------------------------------------------
lvs's avatar
lvs committed
15
/*
16
__kernel __attribute__ ((max_global_work_dim(0)))
lvs's avatar
lvs committed
17
18
*/
__kernel __attribute__ ((reqd_work_group_size(1,1,1)))
19
void Krnl_IntraE(
20
21
 	     __constant     float* restrict KerConstStatic_atom_charges_const,
 	     __constant     char*  restrict KerConstStatic_atom_types_const,
22

23
	     __global const char3* restrict KerConstStatic_intraE_contributors_const,
Leonardo Solis's avatar
Leonardo Solis committed
24
25
26
27
28
29
30

			float 				 DockConst_smooth,
	     __constant     float* restrict KerConstStatic_reqm,
	     __constant     float* restrict KerConstStatic_reqm_hbond,    
	     __constant     uint*  restrict KerConstStatic_atom1_types_reqm,
	     __constant     uint*  restrict KerConstStatic_atom2_types_reqm,  

31
32
33
34
	     __constant     float* restrict KerConstStatic_VWpars_AC_const,
	     __constant     float* restrict KerConstStatic_VWpars_BD_const,
	     __constant     float* restrict KerConstStatic_dspars_S_const,
 	     __constant     float* restrict KerConstStatic_dspars_V_const,
Leonardo Solis's avatar
Leonardo Solis committed
35

36
37
38
39
40
41
			unsigned char                    DockConst_num_of_atoms,
		   	unsigned int                     DockConst_num_of_intraE_contributors,
		  	float                            DockConst_grid_spacing,
			unsigned char                    DockConst_num_of_atypes,
			float                            DockConst_coeff_elec,
			float                            DockConst_qasp,
42
			float                            DockConst_coeff_desolv
43
)
44
{
Leonardo Solis's avatar
Leonardo Solis committed
45
	char active = 0x01;
Leonardo Solis's avatar
Leonardo Solis committed
46

47
48
49
50
51
	__local char3  intraE_contributors_localcache   [MAX_INTRAE_CONTRIBUTORS];
	for (ushort i=0; i<MAX_INTRAE_CONTRIBUTORS; i++) {
		intraE_contributors_localcache [i] = KerConstStatic_intraE_contributors_const [i];	
	}

52
#pragma max_concurrency 32
53
while(active) {
Leonardo Solis's avatar
Leonardo Solis committed
54
55
	char mode;

lvs's avatar
lvs committed
56
/*
57
58
59
60
61
62
63
64
	float3 __attribute__ ((
			      memory,
			      numbanks(2),
			      bankwidth(16),
			      singlepump,
			      numreadports(2),
			      numwriteports(1)
			    )) loc_coords[MAX_NUM_OF_ATOMS];
lvs's avatar
lvs committed
65
66
*/
	float3 loc_coords[MAX_NUM_OF_ATOMS];
67

68
69
70
71
	//printf("BEFORE In INTRA CHANNEL\n");
	// --------------------------------------------------------------
	// Wait for ligand atomic coordinates in channel
	// --------------------------------------------------------------
lvs's avatar
lvs committed
72
/*
Leonardo Solis's avatar
Leonardo Solis committed
73
	char2 actmode = read_channel_altera(chan_Conf2Intrae_actmode);
lvs's avatar
lvs committed
74
75
76
77
*/
	char2 actmode;
	read_pipe_block(chan_Conf2Intrae_actmode, &actmode);
/*
78
	mem_fence(CLK_CHANNEL_MEM_FENCE);
lvs's avatar
lvs committed
79
*/
Leonardo Solis's avatar
Leonardo Solis committed
80

Leonardo Solis's avatar
Leonardo Solis committed
81
82
	active = actmode.x;
	mode   = actmode.y;
83

Leonardo Solis's avatar
Leonardo Solis committed
84
	for (uchar pipe_cnt=0; pipe_cnt<DockConst_num_of_atoms; pipe_cnt+=2) {
lvs's avatar
lvs committed
85
/*
Leonardo Solis's avatar
Leonardo Solis committed
86
		float8 tmp = read_channel_altera(chan_Conf2Intrae_xyz);
lvs's avatar
lvs committed
87
88
89
90
*/
		float8 tmp;
		read_pipe_block(chan_Conf2Intrae_xyz, &tmp);

Leonardo Solis's avatar
Leonardo Solis committed
91
92
93
94
95
96
		float3 tmp1 = {tmp.s0, tmp.s1, tmp.s2};
		float3 tmp2 = {tmp.s4, tmp.s5, tmp.s6};
		loc_coords[pipe_cnt] = tmp1;
		loc_coords[pipe_cnt+1] = tmp2;
	}

97

98
99
	// --------------------------------------------------------------
	//printf("AFTER In INTRA CHANNEL\n");
100

101
102
103
	#if defined (DEBUG_ACTIVE_KERNEL)
	if (active == 0) {printf("	%-20s: %s\n", "Krnl_IntraE", "must be disabled");}
	#endif
Leonardo Solis's avatar
Leonardo Solis committed
104

Leonardo Solis's avatar
Leonardo Solis committed
105
	float intraE = 0.0f;
106

107
108
109
	#if defined (FIXED_POINT_INTRAE)
	// Create shift register to reduce II (initially II=32, unroll-factor=8) 
	// Use fixedpt64 to reduce II=4 (after shift-register) downto II=1
110
111
112
	//float shift_intraE[33];
	fixedpt64 shift_intraE[33];

Leonardo Solis's avatar
Leonardo Solis committed
113
/*
114
	#pragma unroll
Leonardo Solis's avatar
Leonardo Solis committed
115
116
117
*/
	__attribute__((opencl_unroll_hint))
	LOOP_INTRAE_SHIFT_INIT:
118
119
120
121
122
	for (uchar i=0; i<33; i++) {
		//shift_intraE[i] = 0.0f;
		shift_intraE[i] = 0;
	}

123
	#endif
124

125
	// For each intramolecular atom contributor pair
126

127
	//#pragma unroll 10
Leonardo Solis's avatar
Leonardo Solis committed
128
	for (ushort contributor_counter=0; contributor_counter<DockConst_num_of_intraE_contributors; contributor_counter++) {
Leonardo Solis's avatar
Leonardo Solis committed
129

130
131
		char3 ref_intraE_contributors_const;
		ref_intraE_contributors_const = intraE_contributors_localcache[contributor_counter];
Leonardo Solis's avatar
Leonardo Solis committed
132

133
134
		char atom1_id = ref_intraE_contributors_const.x;
		char atom2_id = ref_intraE_contributors_const.y;
Leonardo Solis's avatar
Leonardo Solis committed
135

Leonardo Solis's avatar
Leonardo Solis committed
136
137
138
139
140
		float3 loc_coords_atid1 = loc_coords[atom1_id];
		float3 loc_coords_atid2 = loc_coords[atom2_id];
		float subx = loc_coords_atid1.x - loc_coords_atid2.x;
		float suby = loc_coords_atid1.y - loc_coords_atid2.y;
		float subz = loc_coords_atid1.z - loc_coords_atid2.z;
141

Leonardo Solis's avatar
Leonardo Solis committed
142
143
		//atomic_distance = sqrt(subx*subx + suby*suby + subz*subz)*DockConst_grid_spacing;
		float atomic_distance = sqrt_custom(subx*subx + suby*suby + subz*subz)*DockConst_grid_spacing;
144

Leonardo Solis's avatar
Leonardo Solis committed
145
146
/*
		if (atomic_distance < 1.0f) {
147
			#if defined (DEBUG_KRNL_INTRAE)
Leonardo Solis's avatar
Leonardo Solis committed
148
			printf("\n\nToo low distance (%f) between atoms %u and %u\n", atomic_distance, atom1_id, atom2_id);
149
			#endif
150
			//return HIGHEST_ENERGY;	// Returning maximal value
Leonardo Solis's avatar
Leonardo Solis committed
151
			atomic_distance = 1.0f;
152
		}
Leonardo Solis's avatar
Leonardo Solis committed
153
*/
154

155
		#if defined (DEBUG_KRNL_INTRAE)
156
		printf("\n\nCalculating energy contribution of atoms %u and %u\n", atom1_id+1, atom2_id+1);
Leonardo Solis's avatar
Leonardo Solis committed
157
		printf("Distance: %f\n", atomic_distance);
158
		#endif
Leonardo Solis's avatar
Leonardo Solis committed
159

Leonardo Solis's avatar
Leonardo Solis committed
160
/*
Leonardo Solis's avatar
Leonardo Solis committed
161
162
163
164
		float partialE1;
		float partialE2;
		float partialE3;
		float partialE4;
Leonardo Solis's avatar
Leonardo Solis committed
165
*/
166
		// But only if the distance is less than distance cutoff value and 20.48A (because of the tables)
167
		//if ((distance_leo < 8.0f) && (distance_leo < 20.48f))
Leonardo Solis's avatar
Leonardo Solis committed
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
		float partialE1 = 0.0f;
		float partialE2 = 0.0f;
		float partialE3 = 0.0f;
		float partialE4 = 0.0f;

 		// Getting types ids
		char atom1_typeid = KerConstStatic_atom_types_const [atom1_id];
		char atom2_typeid = KerConstStatic_atom_types_const [atom2_id];

		// Getting optimum pair distance (opt_distance) from reqm and reqm_hbond
		// reqm: equilibrium internuclear separation 
		//       (sum of the vdW radii of two like atoms (A)) in the case of vdW
		// reqm_hbond: equilibrium internuclear separation
		// 	 (sum of the vdW radii of two like atoms (A)) in the case of hbond 
		float opt_distance;

		uint atom1_type_vdw_hb = KerConstStatic_atom1_types_reqm [atom1_typeid];
     	        uint atom2_type_vdw_hb = KerConstStatic_atom2_types_reqm [atom2_typeid];

		if (ref_intraE_contributors_const.z == 1)	// H-bond
188
		{
Leonardo Solis's avatar
Leonardo Solis committed
189
190
191
192
193
194
			opt_distance = KerConstStatic_reqm_hbond [atom1_type_vdw_hb] + KerConstStatic_reqm_hbond [atom2_type_vdw_hb];
		}
		else	// Van der Waals
		{
			opt_distance = 0.5f*(KerConstStatic_reqm [atom1_type_vdw_hb] + KerConstStatic_reqm [atom2_type_vdw_hb]);
		}
195

Leonardo Solis's avatar
Leonardo Solis committed
196
197
198
199
		// Getting smoothed distance
		// smoothed_distance = function(atomic_distance, opt_distance)
		float smoothed_distance;
		float delta_distance = 0.5f*DockConst_smooth;
200

Leonardo Solis's avatar
Leonardo Solis committed
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
		/*printf("delta_distance: %f\n", delta_distance);*/

		if (atomic_distance <= (opt_distance - delta_distance)) {
			smoothed_distance = atomic_distance + delta_distance;
		}
		else if (atomic_distance < (opt_distance + delta_distance)) {
			smoothed_distance = opt_distance;
		}
		else { // else if (atomic_distance >= (opt_distance + delta_distance))
			smoothed_distance = atomic_distance - delta_distance;
		}

		float distance_pow_2  = atomic_distance*atomic_distance; 

		float smoothed_distance_pow_2 = smoothed_distance*smoothed_distance; 
		float inverse_smoothed_distance_pow_2  = native_divide(1.0f, smoothed_distance_pow_2);
		float inverse_smoothed_distance_pow_4  = inverse_smoothed_distance_pow_2 * inverse_smoothed_distance_pow_2;
		float inverse_smoothed_distance_pow_6  = inverse_smoothed_distance_pow_4 * inverse_smoothed_distance_pow_2;
		float inverse_smoothed_distance_pow_10 = inverse_smoothed_distance_pow_6 * inverse_smoothed_distance_pow_4;
		float inverse_smoothed_distance_pow_12 = inverse_smoothed_distance_pow_6 * inverse_smoothed_distance_pow_6;

		if (atomic_distance < 8.0f) 
		{
224
			// Calculating van der Waals / hydrogen bond term
Leonardo Solis's avatar
Leonardo Solis committed
225
			partialE1 = KerConstStatic_VWpars_AC_const [atom1_typeid*DockConst_num_of_atypes+atom2_typeid]*inverse_smoothed_distance_pow_12;
226

227
228
			float tmp_pE2 = KerConstStatic_VWpars_BD_const [atom1_typeid*DockConst_num_of_atypes+atom2_typeid];

229
			if (ref_intraE_contributors_const.z == 1)	// H-bond
Leonardo Solis's avatar
Leonardo Solis committed
230
				partialE2 = tmp_pE2 * inverse_smoothed_distance_pow_10;
231
			else	// Van der Waals
Leonardo Solis's avatar
Leonardo Solis committed
232
233
				partialE2 = tmp_pE2 * inverse_smoothed_distance_pow_6;
		} // End of if: if (dist < dcutoff)	
234

Leonardo Solis's avatar
Leonardo Solis committed
235
236
		// Calculating electrostatic term
		partialE3 = native_divide(  (DockConst_coeff_elec*KerConstStatic_atom_charges_const[atom1_id]*KerConstStatic_atom_charges_const[atom2_id]) , (atomic_distance*(-8.5525f + native_divide(86.9525f, (1.0f + 7.7839f*native_exp(-0.3154f*atomic_distance)))))       );
237

Leonardo Solis's avatar
Leonardo Solis committed
238
239
240
241
242
		// Calculating desolvation term
		partialE4 = (
			  (KerConstStatic_dspars_S_const [atom1_typeid] + DockConst_qasp*fabs(KerConstStatic_atom_charges_const[atom1_id])) * KerConstStatic_dspars_V_const [atom2_typeid] + 
			  (KerConstStatic_dspars_S_const [atom2_typeid] + DockConst_qasp*fabs(KerConstStatic_atom_charges_const[atom2_id])) * KerConstStatic_dspars_V_const [atom1_typeid]) * 
			 DockConst_coeff_desolv*native_exp(-0.0386f*distance_pow_2);
Leonardo Solis's avatar
Leonardo Solis committed
243
	
244
		#if defined (FIXED_POINT_INTRAE)
245
246
247
248
249
250
		//shift_intraE[32] = shift_intraE[0] + partialE1 + partialE2 + partialE3 + partialE4;
		shift_intraE[32] = shift_intraE[0] + fixedpt64_fromfloat(partialE1) + 
						     fixedpt64_fromfloat(partialE2) + 
						     fixedpt64_fromfloat(partialE3) + 
						     fixedpt64_fromfloat(partialE4);

Leonardo Solis's avatar
Leonardo Solis committed
251
/*
252
		#pragma unroll
Leonardo Solis's avatar
Leonardo Solis committed
253
254
255
*/
		__attribute__((opencl_unroll_hint))
		LOOP_INTRAE_SHIFT:
256
257
258
		for (uchar j=0; j<32; j++) {
			shift_intraE[j] = shift_intraE[j+1];
		}
259
		#else
Leonardo Solis's avatar
Leonardo Solis committed
260
		intraE += partialE1 + partialE2 + partialE3 + partialE4;
261
		#endif
Leonardo Solis's avatar
Leonardo Solis committed
262
	
Leonardo Solis's avatar
Leonardo Solis committed
263
	} // End of contributor_counter for-loop
264

265
	#if defined (FIXED_POINT_INTRAE)
266
267
	fixedpt64 fixpt_intraE = 0;

Leonardo Solis's avatar
Leonardo Solis committed
268
/*
269
	#pragma unroll
Leonardo Solis's avatar
Leonardo Solis committed
270
271
272
*/
	__attribute__((opencl_unroll_hint))
	LOOP_INTRAE_SHIFT_RED:
273
274
275
276
277
	for (uchar j=0; j<32; j++) {
		//intraE += shift_intraE[j];
		fixpt_intraE += shift_intraE[j];
	}
	intraE = fixedpt64_tofloat(fixpt_intraE);
278
	#endif
279

280
281
282
	// --------------------------------------------------------------
	// Send intramolecular energy to channel
	// --------------------------------------------------------------
283
	switch (mode) {
lvs's avatar
lvs committed
284
/*
285
286
		case 'I':  write_channel_altera(chan_Intrae2StoreIC_intrae, intraE);     break;
		case 'G':  write_channel_altera(chan_Intrae2StoreGG_intrae, intraE);     break;
287
288
289
		case 0x01: write_channel_altera(chan_Intrae2StoreLS_LS1_intrae, intraE); break;
		case 0x02: write_channel_altera(chan_Intrae2StoreLS_LS2_intrae, intraE); break;
		case 0x03: write_channel_altera(chan_Intrae2StoreLS_LS3_intrae, intraE); break;
Leonardo Solis's avatar
Leonardo Solis committed
290
291
		case 0x04: write_channel_altera(chan_Intrae2StoreLS_LS4_intrae, intraE); break;
		case 0x05: write_channel_altera(chan_Intrae2StoreLS_LS5_intrae, intraE); break;
292
293
294
295
		case 0x06: write_channel_altera(chan_Intrae2StoreLS_LS6_intrae, intraE); break;
		case 0x07: write_channel_altera(chan_Intrae2StoreLS_LS7_intrae, intraE); break;
		case 0x08: write_channel_altera(chan_Intrae2StoreLS_LS8_intrae, intraE); break;
		case 0x09: write_channel_altera(chan_Intrae2StoreLS_LS9_intrae, intraE); break;
lvs's avatar
lvs committed
296
297
298
299
300
301
302
303
304
305
306
307
*/
		case 'I':  write_pipe_block(chan_Intrae2StoreIC_intrae, &intraE);     break;
		case 'G':  write_pipe_block(chan_Intrae2StoreGG_intrae, &intraE);     break;
		case 0x01: write_pipe_block(chan_Intrae2StoreLS_LS1_intrae, &intraE); break;
		case 0x02: write_pipe_block(chan_Intrae2StoreLS_LS2_intrae, &intraE); break;
		case 0x03: write_pipe_block(chan_Intrae2StoreLS_LS3_intrae, &intraE); break;
		case 0x04: write_pipe_block(chan_Intrae2StoreLS_LS4_intrae, &intraE); break;
		case 0x05: write_pipe_block(chan_Intrae2StoreLS_LS5_intrae, &intraE); break;
		case 0x06: write_pipe_block(chan_Intrae2StoreLS_LS6_intrae, &intraE); break;
		case 0x07: write_pipe_block(chan_Intrae2StoreLS_LS7_intrae, &intraE); break;
		case 0x08: write_pipe_block(chan_Intrae2StoreLS_LS8_intrae, &intraE); break;
		case 0x09: write_pipe_block(chan_Intrae2StoreLS_LS9_intrae, &intraE); break;
308
	}
309
	// --------------------------------------------------------------
310

Leonardo Solis's avatar
Leonardo Solis committed
311
312
313
314
	#if defined (DEBUG_KRNL_INTRAE)
	printf("AFTER Out INTRAE CHANNEL\n");
	#endif

315
} // End of while(active)
316
317
318
319

	#if defined (DEBUG_ACTIVE_KERNEL)
	printf("	%-20s: %s\n", "Krnl_IntraE", "disabled");
	#endif
320
321
322
}
// --------------------------------------------------------------------------
// --------------------------------------------------------------------------