Krnl_IntraE.cl 9.82 KB
Newer Older
1
// sqrt7 ////https://www.codeproject.com/Articles/69941/Best-Square-Root-Method-Algorithm-Function-Precisi
2
__attribute__((always_inline))
3
4
5
6
7
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 	
8
	return as_float(i);	// return *(float*) &i; 
9
10
}  

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

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

			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,  

29
30
31
32
	     __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
33

34
35
36
37
38
39
			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,
40
			float                            DockConst_coeff_desolv
41
)
42
{
Leonardo Solis's avatar
Leonardo Solis committed
43
	char active = 0x01;
Leonardo Solis's avatar
Leonardo Solis committed
44

45
	__local char3  intraE_contributors_localcache   [MAX_INTRAE_CONTRIBUTORS];
46
47

	__attribute__((xcl_pipeline_loop))
Leonardo Solis's avatar
Leonardo Solis committed
48
	LOOP_FOR_INTRAE_CONTRIBUTORS:
49
50
51
52
	for (ushort i=0; i<MAX_INTRAE_CONTRIBUTORS; i++) {
		intraE_contributors_localcache [i] = KerConstStatic_intraE_contributors_const [i];	
	}

53
__attribute__((xcl_pipeline_loop))
Leonardo Solis's avatar
Leonardo Solis committed
54
LOOP_WHILE_INTRAE_MAIN:
55
while(active) {
Leonardo Solis's avatar
Leonardo Solis committed
56
57
	char mode;

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

70
71
72
73
	//printf("BEFORE In INTRA CHANNEL\n");
	// --------------------------------------------------------------
	// Wait for ligand atomic coordinates in channel
	// --------------------------------------------------------------
Leonardo Solis's avatar
Leonardo Solis committed
74

75
	char actmode;
lvs's avatar
lvs committed
76
77
	read_pipe_block(chan_Conf2Intrae_actmode, &actmode);
/*
78
	mem_fence(CLK_CHANNEL_MEM_FENCE);
lvs's avatar
lvs committed
79
*/
80
81
	active = actmode;
	mode   = actmode;
82

83
	__attribute__((xcl_pipeline_loop))
Leonardo Solis's avatar
Leonardo Solis committed
84
	LOOP_FOR_INTRAE_READ_XYZ:
Leonardo Solis's avatar
Leonardo Solis committed
85
	for (uchar pipe_cnt=0; pipe_cnt<DockConst_num_of_atoms; pipe_cnt+=2) {
lvs's avatar
lvs committed
86
87
88
		float8 tmp;
		read_pipe_block(chan_Conf2Intrae_xyz, &tmp);

Leonardo Solis's avatar
Leonardo Solis committed
89
90
91
92
93
94
		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;
	}

95
96
	// --------------------------------------------------------------
	//printf("AFTER In INTRA CHANNEL\n");
97

98
99
100
	#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
101

Leonardo Solis's avatar
Leonardo Solis committed
102
	float intraE = 0.0f;
103

104

105
	// For each intramolecular atom contributor pair
106

107
	__attribute__((xcl_pipeline_loop))
Leonardo Solis's avatar
Leonardo Solis committed
108
	LOOP_FOR_INTRAE_MAIN:
Leonardo Solis's avatar
Leonardo Solis committed
109
	for (ushort contributor_counter=0; contributor_counter<DockConst_num_of_intraE_contributors; contributor_counter++) {
Leonardo Solis's avatar
Leonardo Solis committed
110

111
112
		char3 ref_intraE_contributors_const;
		ref_intraE_contributors_const = intraE_contributors_localcache[contributor_counter];
Leonardo Solis's avatar
Leonardo Solis committed
113

114
115
		char atom1_id = ref_intraE_contributors_const.x;
		char atom2_id = ref_intraE_contributors_const.y;
Leonardo Solis's avatar
Leonardo Solis committed
116

Leonardo Solis's avatar
Leonardo Solis committed
117
118
119
120
121
		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;
122

Leonardo Solis's avatar
Leonardo Solis committed
123
124
		//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;
125

Leonardo Solis's avatar
Leonardo Solis committed
126
127
/*
		if (atomic_distance < 1.0f) {
128
			#if defined (DEBUG_KRNL_INTRAE)
Leonardo Solis's avatar
Leonardo Solis committed
129
			printf("\n\nToo low distance (%f) between atoms %u and %u\n", atomic_distance, atom1_id, atom2_id);
130
			#endif
131
			//return HIGHEST_ENERGY;	// Returning maximal value
Leonardo Solis's avatar
Leonardo Solis committed
132
			atomic_distance = 1.0f;
133
		}
Leonardo Solis's avatar
Leonardo Solis committed
134
*/
135

136
		#if defined (DEBUG_KRNL_INTRAE)
137
		printf("\n\nCalculating energy contribution of atoms %u and %u\n", atom1_id+1, atom2_id+1);
Leonardo Solis's avatar
Leonardo Solis committed
138
		printf("Distance: %f\n", atomic_distance);
139
		#endif
Leonardo Solis's avatar
Leonardo Solis committed
140

Leonardo Solis's avatar
Leonardo Solis committed
141
/*
Leonardo Solis's avatar
Leonardo Solis committed
142
143
144
145
		float partialE1;
		float partialE2;
		float partialE3;
		float partialE4;
Leonardo Solis's avatar
Leonardo Solis committed
146
*/
147
		// But only if the distance is less than distance cutoff value and 20.48A (because of the tables)
148
		//if ((distance_leo < 8.0f) && (distance_leo < 20.48f))
Leonardo Solis's avatar
Leonardo Solis committed
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
		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
169
		{
Leonardo Solis's avatar
Leonardo Solis committed
170
171
172
173
174
175
			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]);
		}
176

Leonardo Solis's avatar
Leonardo Solis committed
177
178
179
180
		// Getting smoothed distance
		// smoothed_distance = function(atomic_distance, opt_distance)
		float smoothed_distance;
		float delta_distance = 0.5f*DockConst_smooth;
181

Leonardo Solis's avatar
Leonardo Solis committed
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
		/*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;

203
204
		// Calculating energy contributions
		// Cuttoff1: internuclear-distance at 8A only for vdw and hbond.
Leonardo Solis's avatar
Leonardo Solis committed
205
206
		if (atomic_distance < 8.0f) 
		{
207
			// Calculating van der Waals / hydrogen bond term
208
			partialE1 += KerConstStatic_VWpars_AC_const [atom1_typeid*DockConst_num_of_atypes+atom2_typeid]*inverse_smoothed_distance_pow_12;
209

210
211
			float tmp_pE2 = KerConstStatic_VWpars_BD_const [atom1_typeid*DockConst_num_of_atypes+atom2_typeid];

212
			if (ref_intraE_contributors_const.z == 1)	// H-bond
213
				partialE2 -= tmp_pE2 * inverse_smoothed_distance_pow_10;
214
			else	// Van der Waals
215
216
				partialE2 -= tmp_pE2 * inverse_smoothed_distance_pow_6;
		} // if cuttoff1 - internuclear-distance at 8A
217

218
219
220
221
222
223
224
225
226
227
228
229
230
		// Calculating energy contributions
		// Cuttoff2: internuclear-distance at 20.48A only for el and sol.
		if (atomic_distance < 20.48f)
		{
			// 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)))))       );

			// 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);
		} // if cuttoff2 - internuclear-distance at 20.48A
Leonardo Solis's avatar
Leonardo Solis committed
231
232
233
	
		intraE += partialE1 + partialE2 + partialE3 + partialE4;
	
Leonardo Solis's avatar
Leonardo Solis committed
234
	} // End of contributor_counter for-loop
235

236
237
238
	// --------------------------------------------------------------
	// Send intramolecular energy to channel
	// --------------------------------------------------------------
239
	switch (mode) {
Leonardo Solis's avatar
Leonardo Solis committed
240
241
		case 'I':  write_pipe_block(chan_Intrae2StoreIC_intrae,     &intraE); break;
		case 'G':  write_pipe_block(chan_Intrae2StoreGG_intrae,     &intraE); break;
242
	}
243
	// --------------------------------------------------------------
244

Leonardo Solis's avatar
Leonardo Solis committed
245
246
247
248
	#if defined (DEBUG_KRNL_INTRAE)
	printf("AFTER Out INTRAE CHANNEL\n");
	#endif

249
} // End of while(active)
250
251
252
253

	#if defined (DEBUG_ACTIVE_KERNEL)
	printf("	%-20s: %s\n", "Krnl_IntraE", "disabled");
	#endif
254
255
256
}
// --------------------------------------------------------------------------
// --------------------------------------------------------------------------