Tinker9 70bd052 (Thu Nov 9 12:11:35 2023 -0800)
Loading...
Searching...
No Matches
geom.h
1#pragma once
2#include "ff/image.h"
3#include "math/const.h"
4#include "math/libfunc.h"
5#include "seq/add.h"
6#include "seq/seq.h"
7
8namespace tinker {
9#pragma acc routine seq
10template <class Ver>
13 real& restrict vzx, real& restrict vyy, real& restrict vzy,
14 real& restrict vzz,
15
17
18 int i, const int (*restrict igfix)[2], const real (*restrict gfix)[3],
19
20 const real* restrict x, const real* restrict y, const real* restrict z,
21 const double* restrict mass, const int* restrict molec,
22 const int (*restrict igrp)[2], const int* restrict kgrp,
23 const double* restrict grpmass, TINKER_IMAGE_PARAMS)
24{
25 constexpr bool do_e = Ver::e;
26 constexpr bool do_g = Ver::g;
27 constexpr bool do_v = Ver::v;
28
29 int ia = igfix[i][0];
30 int ib = igfix[i][1];
31 int ja1 = igrp[ia][0];
32 int ja2 = igrp[ia][1];
33 int jb1 = igrp[ib][0];
34 int jb2 = igrp[ib][1];
35
36 // image() calls are necessary here because atoms in the same group
37 // do not necessarily belong to the same molecule
38
39 real weigha = REAL_MAX((real)1, grpmass[ia]);
40 int ka1 = kgrp[ja1];
41 real xa1 = x[ka1], ya1 = y[ka1], za1 = z[ka1];
42 real xacm = xa1 * weigha, yacm = ya1 * weigha, zacm = za1 * weigha;
43 #pragma acc loop seq
44 for (int j = ja1 + 1; j < ja2; ++j) {
45 int k = kgrp[j];
46 real weigh = mass[k];
47 real xr = x[k] - xa1;
48 real yr = y[k] - ya1;
49 real zr = z[k] - za1;
50 if (molec[ka1] != molec[k]) image(xr, yr, zr);
51 xacm += xr * weigh;
52 yacm += yr * weigh;
53 zacm += zr * weigh;
54 }
55 weigha = REAL_RECIP(weigha);
56
57 real weighb = REAL_MAX((real)1, grpmass[ib]);
58 int kb1 = kgrp[jb1];
59 real xb1 = x[kb1], yb1 = y[kb1], zb1 = z[kb1];
60 real xbcm = xb1 * weighb, ybcm = yb1 * weighb, zbcm = zb1 * weighb;
61 #pragma acc loop seq
62 for (int j = jb1 + 1; j < jb2; ++j) {
63 int k = kgrp[j];
64 real weigh = mass[k];
65 real xr = x[k] - xb1;
66 real yr = y[k] - yb1;
67 real zr = z[k] - zb1;
68 if (molec[kb1] != molec[k]) image(xr, yr, zr);
69 xbcm += xr * weigh;
70 ybcm += yr * weigh;
71 zbcm += zr * weigh;
72 }
73 weighb = REAL_RECIP(weighb);
74
75 real xr = xacm * weigha - xbcm * weighb;
76 real yr = yacm * weigha - ybcm * weighb;
77 real zr = zacm * weigha - zbcm * weighb;
78
79 bool intermol = molec[kgrp[ja1]] != molec[kgrp[jb1]];
80 if (intermol) image(xr, yr, zr);
81
82 real r = REAL_SQRT(xr * xr + yr * yr + zr * zr);
83 real force = gfix[i][0];
84 real gf1 = gfix[i][1];
85 real gf2 = gfix[i][2];
86 real target = (r < gf1 ? gf1 : (r > gf2 ? gf2 : r));
87 real dt = r - target;
88
89 if CONSTEXPR (do_e) {
90 real dt2 = dt * dt;
91 e = force * dt2;
92 }
93 if CONSTEXPR (do_g) {
94 real rinv = (r == 0 ? 1 : REAL_RECIP(r));
95 real de = 2 * force * dt * rinv;
96 real dedx = de * xr;
97 real dedy = de * yr;
98 real dedz = de * zr;
99
100 #pragma acc loop seq
101 for (int j = ja1; j < ja2; ++j) {
102 int k = kgrp[j];
103 real ratio = mass[k] * weigha;
104 atomic_add(dedx * ratio, degx, k);
105 atomic_add(dedy * ratio, degy, k);
106 atomic_add(dedz * ratio, degz, k);
107 }
108 #pragma acc loop seq
109 for (int j = jb1; j < jb2; ++j) {
110 int k = kgrp[j];
111 real ratio = mass[k] * weighb;
112 atomic_add(-dedx * ratio, degx, k);
113 atomic_add(-dedy * ratio, degy, k);
114 atomic_add(-dedz * ratio, degz, k);
115 }
116 if CONSTEXPR (do_v) {
117 vxx = xr * dedx;
118 vyx = yr * dedx;
119 vzx = zr * dedx;
120 vyy = yr * dedy;
121 vzy = zr * dedy;
122 vzz = zr * dedz;
123 }
124 }
125}
126
127#pragma acc routine seq
128template <class Ver>
131 real& restrict vzx, real& restrict vyy, real& restrict vzy,
132 real& restrict vzz,
133
135
136 int i, const int (*restrict idfix)[2], const real (*restrict dfix)[3],
137
138 const real* restrict x, const real* restrict y, const real* restrict z,
139 const int* restrict molec, TINKER_IMAGE_PARAMS)
140{
141 constexpr bool do_e = Ver::e;
142 constexpr bool do_g = Ver::g;
143 constexpr bool do_v = Ver::v;
144
145 int ia = idfix[i][0];
146 int ib = idfix[i][1];
147 real force = dfix[i][0];
148 real df1 = dfix[i][1];
149 real df2 = dfix[i][2];
150
151 real xr = x[ia] - x[ib];
152 real yr = y[ia] - y[ib];
153 real zr = z[ia] - z[ib];
154 bool intermol = molec[ia] != molec[ib];
155 if (intermol) image(xr, yr, zr);
156 real r = REAL_SQRT(xr * xr + yr * yr + zr * zr);
157 real target = r;
158 if (r < df1) target = df1;
159 if (r > df2) target = df2;
160 real dt = r - target;
161
162 if CONSTEXPR (do_e) {
163 real dt2 = dt * dt;
164 e = force * dt2;
165 }
166 if CONSTEXPR (do_g) {
167 real rinv = (r == 0 ? 1 : REAL_RECIP(r));
168 real de = 2 * force * dt * rinv;
169 real dedx = de * xr;
170 real dedy = de * yr;
171 real dedz = de * zr;
172 atomic_add(dedx, degx, ia);
173 atomic_add(dedy, degy, ia);
174 atomic_add(dedz, degz, ia);
175 atomic_add(-dedx, degx, ib);
176 atomic_add(-dedy, degy, ib);
177 atomic_add(-dedz, degz, ib);
178 if CONSTEXPR (do_v) {
179 vxx = xr * dedx;
180 vyx = yr * dedx;
181 vzx = zr * dedx;
182 vyy = yr * dedy;
183 vzy = zr * dedy;
184 vzz = zr * dedz;
185 }
186 }
187}
188
189#pragma acc routine seq
190template <class Ver>
193 real& restrict vzx, real& restrict vyy, real& restrict vzy,
194 real& restrict vzz,
195
197
198 int i, const int (*restrict iafix)[3], const real (*restrict afix)[3],
199
200 const real* restrict x, const real* restrict y, const real* restrict z)
201{
202 constexpr bool do_e = Ver::e;
203 constexpr bool do_g = Ver::g;
204 constexpr bool do_v = Ver::v;
205 if CONSTEXPR (do_e) e = 0;
206 if CONSTEXPR (do_v) vxx = 0, vyx = 0, vzx = 0, vyy = 0, vzy = 0, vzz = 0;
207
208 int ia = iafix[i][0];
209 int ib = iafix[i][1];
210 int ic = iafix[i][2];
211 real force = afix[i][0];
212 real af1 = afix[i][1];
213 real af2 = afix[i][2];
214
215 real xia = x[ia];
216 real yia = y[ia];
217 real zia = z[ia];
218 real xib = x[ib];
219 real yib = y[ib];
220 real zib = z[ib];
221 real xic = x[ic];
222 real yic = y[ic];
223 real zic = z[ic];
224
225 real xab = xia - xib;
226 real yab = yia - yib;
227 real zab = zia - zib;
228 real xcb = xic - xib;
229 real ycb = yic - yib;
230 real zcb = zic - zib;
231 real rab2 = xab * xab + yab * yab + zab * zab;
232 rab2 = REAL_MAX(rab2, (real)0.0001);
233 real rcb2 = xcb * xcb + ycb * ycb + zcb * zcb;
234 rcb2 = REAL_MAX(rcb2, (real)0.0001);
235
236 real xp = ycb * zab - zcb * yab;
237 real yp = zcb * xab - xcb * zab;
238 real zp = xcb * yab - ycb * xab;
239 real rp = REAL_SQRT(xp * xp + yp * yp + zp * zp);
240 rp = REAL_MAX(rp, (real)0.0001);
241 real dot = xab * xcb + yab * ycb + zab * zcb;
242 real cosine = dot * REAL_RSQRT(rab2 * rcb2);
243 cosine = REAL_MIN((real)1, REAL_MAX((real)-1, cosine));
244 real angle = radian * REAL_ACOS(cosine);
245
246 real target = angle;
247 if (angle < af1) target = af1;
248 if (angle > af2) target = af2;
249 real dt = (angle - target) * _1radian;
250
251 if CONSTEXPR (do_e) {
252 real dt2 = dt * dt;
253 e = force * dt2;
254 }
255 if CONSTEXPR (do_g) {
256 real deddt = 2 * force * dt;
257 real terma = -deddt * REAL_RECIP(rab2 * rp);
258 real termc = deddt * REAL_RECIP(rcb2 * rp);
259 real dedxia = terma * (yab * zp - zab * yp);
260 real dedyia = terma * (zab * xp - xab * zp);
261 real dedzia = terma * (xab * yp - yab * xp);
262 real dedxic = termc * (ycb * zp - zcb * yp);
263 real dedyic = termc * (zcb * xp - xcb * zp);
264 real dedzic = termc * (xcb * yp - ycb * xp);
265 real dedxib = -dedxia - dedxic;
266 real dedyib = -dedyia - dedyic;
267 real dedzib = -dedzia - dedzic;
268
269 atomic_add(dedxia, degx, ia);
270 atomic_add(dedyia, degy, ia);
271 atomic_add(dedzia, degz, ia);
272 atomic_add(dedxib, degx, ib);
273 atomic_add(dedyib, degy, ib);
274 atomic_add(dedzib, degz, ib);
275 atomic_add(dedxic, degx, ic);
276 atomic_add(dedyic, degy, ic);
277 atomic_add(dedzic, degz, ic);
278 if CONSTEXPR (do_v) {
279 vxx = xab * dedxia + xcb * dedxic;
280 vyx = yab * dedxia + ycb * dedxic;
281 vzx = zab * dedxia + zcb * dedxic;
282 vyy = yab * dedyia + ycb * dedyic;
283 vzy = zab * dedyia + zcb * dedyic;
284 vzz = zab * dedzia + zcb * dedzic;
285 }
286 }
287}
288
289#pragma acc routine seq
290template <class Ver>
293 real& restrict vzx, real& restrict vyy, real& restrict vzy,
294 real& restrict vzz,
295
297
298 int i, const int (*restrict itfix)[4], const real (*restrict tfix)[3],
299
300 const real* restrict x, const real* restrict y, const real* restrict z)
301{
302 constexpr bool do_e = Ver::e;
303 constexpr bool do_g = Ver::g;
304 constexpr bool do_v = Ver::v;
305 if CONSTEXPR (do_e) e = 0;
306 if CONSTEXPR (do_v) vxx = 0, vyx = 0, vzx = 0, vyy = 0, vzy = 0, vzz = 0;
307
308 int ia = itfix[i][0];
309 int ib = itfix[i][1];
310 int ic = itfix[i][2];
311 int id = itfix[i][3];
312 real force = tfix[i][0];
313 real tf1 = tfix[i][1];
314 real tf2 = tfix[i][2];
315
316 real xia = x[ia];
317 real yia = y[ia];
318 real zia = z[ia];
319 real xib = x[ib];
320 real yib = y[ib];
321 real zib = z[ib];
322 real xic = x[ic];
323 real yic = y[ic];
324 real zic = z[ic];
325 real xid = x[id];
326 real yid = y[id];
327 real zid = z[id];
328 real xba = xib - xia;
329 real yba = yib - yia;
330 real zba = zib - zia;
331 real xcb = xic - xib;
332 real ycb = yic - yib;
333 real zcb = zic - zib;
334 real rcb = REAL_SQRT(xcb * xcb + ycb * ycb + zcb * zcb);
335 rcb = REAL_MAX(rcb, (real)0.0001);
336 real xdc = xid - xic;
337 real ydc = yid - yic;
338 real zdc = zid - zic;
339
340 real xt = yba * zcb - ycb * zba;
341 real yt = zba * xcb - zcb * xba;
342 real zt = xba * ycb - xcb * yba;
343 real xu = ycb * zdc - ydc * zcb;
344 real yu = zcb * xdc - zdc * xcb;
345 real zu = xcb * ydc - xdc * ycb;
346 real xtu = yt * zu - yu * zt;
347 real ytu = zt * xu - zu * xt;
348 real ztu = xt * yu - xu * yt;
349 real rt2 = xt * xt + yt * yt + zt * zt;
350 rt2 = REAL_MAX(rt2, (real)0.0001);
351 real ru2 = xu * xu + yu * yu + zu * zu;
352 ru2 = REAL_MAX(ru2, (real)0.0001);
353 real rtru = REAL_SQRT(rt2 * ru2);
354
355 real cosine = (xt * xu + yt * yu + zt * zu) * REAL_RECIP(rtru);
356 real sine = (xcb * xtu + ycb * ytu + zcb * ztu) * REAL_RECIP(rcb * rtru);
357 cosine = REAL_MIN((real)1.0, REAL_MAX((real)-1.0, cosine));
358 real angle = radian * REAL_ACOS(cosine);
359 if (sine < 0) angle = -angle;
360 real target;
361 if (angle > tf1 and angle < tf2)
362 target = angle;
363 else if (angle > tf1 and tf1 > tf2)
364 target = angle;
365 else if (angle < tf2 and tf1 > tf2)
366 target = angle;
367 else {
368 real t1 = angle - tf1;
369 real t2 = angle - tf2;
370 if (t1 > 180) t1 -= 360;
371 if (t1 < -180) t1 += 360;
372 if (t2 > 180) t2 -= 360;
373 if (t2 < -180) t2 += 360;
374 if (REAL_ABS(t1) < REAL_ABS(t2))
375 target = tf1;
376 else
377 target = tf2;
378 }
379 real dt = angle - target;
380 if (dt > 180) dt -= 360;
381 if (dt < -180) dt += 360;
382
383 dt *= _1radian;
384 if CONSTEXPR (do_e) {
385 real dt2 = dt * dt;
386 e = force * dt2;
387 }
388 if CONSTEXPR (do_g) {
389 real dedphi = 2 * force * dt;
390
391 real xca = xic - xia;
392 real yca = yic - yia;
393 real zca = zic - zia;
394 real xdb = xid - xib;
395 real ydb = yid - yib;
396 real zdb = zid - zib;
397
398 real rt_inv = REAL_RECIP(rt2 * rcb);
399 real ru_inv = REAL_RECIP(ru2 * rcb);
400 real dedxt = dedphi * (yt * zcb - ycb * zt) * rt_inv;
401 real dedyt = dedphi * (zt * xcb - zcb * xt) * rt_inv;
402 real dedzt = dedphi * (xt * ycb - xcb * yt) * rt_inv;
403 real dedxu = -dedphi * (yu * zcb - ycb * zu) * ru_inv;
404 real dedyu = -dedphi * (zu * xcb - zcb * xu) * ru_inv;
405 real dedzu = -dedphi * (xu * ycb - xcb * yu) * ru_inv;
406
407 real dedxia = zcb * dedyt - ycb * dedzt;
408 real dedyia = xcb * dedzt - zcb * dedxt;
409 real dedzia = ycb * dedxt - xcb * dedyt;
410 real dedxib = yca * dedzt - zca * dedyt + zdc * dedyu - ydc * dedzu;
411 real dedyib = zca * dedxt - xca * dedzt + xdc * dedzu - zdc * dedxu;
412 real dedzib = xca * dedyt - yca * dedxt + ydc * dedxu - xdc * dedyu;
413 real dedxic = zba * dedyt - yba * dedzt + ydb * dedzu - zdb * dedyu;
414 real dedyic = xba * dedzt - zba * dedxt + zdb * dedxu - xdb * dedzu;
415 real dedzic = yba * dedxt - xba * dedyt + xdb * dedyu - ydb * dedxu;
416 real dedxid = zcb * dedyu - ycb * dedzu;
417 real dedyid = xcb * dedzu - zcb * dedxu;
418 real dedzid = ycb * dedxu - xcb * dedyu;
419
420 atomic_add(dedxia, degx, ia);
421 atomic_add(dedyia, degy, ia);
422 atomic_add(dedzia, degz, ia);
423 atomic_add(dedxib, degx, ib);
424 atomic_add(dedyib, degy, ib);
425 atomic_add(dedzib, degz, ib);
426 atomic_add(dedxic, degx, ic);
427 atomic_add(dedyic, degy, ic);
428 atomic_add(dedzic, degz, ic);
429 atomic_add(dedxid, degx, id);
430 atomic_add(dedyid, degy, id);
431 atomic_add(dedzid, degz, id);
432
433 if CONSTEXPR (do_v) {
434 vxx = xcb * (dedxic + dedxid) - xba * dedxia + xdc * dedxid;
435 vyx = ycb * (dedxic + dedxid) - yba * dedxia + ydc * dedxid;
436 vzx = zcb * (dedxic + dedxid) - zba * dedxia + zdc * dedxid;
437 vyy = ycb * (dedyic + dedyid) - yba * dedyia + ydc * dedyid;
438 vzy = zcb * (dedyic + dedyid) - zba * dedyia + zdc * dedyid;
439 vzz = zcb * (dedzic + dedzid) - zba * dedzia + zdc * dedzid;
440 }
441 }
442}
443
444#pragma acc routine seq
445template <class Ver>
448 real& restrict vzx, real& restrict vyy, real& restrict vzy,
449 real& restrict vzz,
450
452
453 int i, const int* restrict ipfix, const int (*restrict kpfix)[3],
454 const real* restrict xpfix, const real* restrict ypfix,
455 const real* restrict zpfix, const real (*restrict pfix)[2],
456
457 const real* restrict x, const real* restrict y, const real* restrict z,
459{
460 constexpr bool do_e = Ver::e;
461 constexpr bool do_g = Ver::g;
462 constexpr bool do_v = Ver::v;
463
464 int ia = ipfix[i];
465 real force = pfix[i][0];
466 real radius = pfix[i][1];
467 real xr = 0, yr = 0, zr = 0;
468 if (kpfix[i][0]) xr = x[ia] - xpfix[i];
469 if (kpfix[i][1]) yr = y[ia] - ypfix[i];
470 if (kpfix[i][2]) zr = z[ia] - zpfix[i];
471 image(xr, yr, zr);
472 real r = REAL_SQRT(xr * xr + yr * yr + zr * zr);
473 real dt = REAL_MAX((real)0, r - radius);
474
475 if CONSTEXPR (do_e) {
476 real dt2 = dt * dt;
477 e = force * dt2;
478 }
479 if CONSTEXPR (do_g) {
480 real rinv = (r == 0 ? 1 : REAL_RECIP(r));
481 real de = 2 * force * dt * rinv;
482 real dedx = de * xr;
483 real dedy = de * yr;
484 real dedz = de * zr;
485 atomic_add(dedx, degx, ia);
486 atomic_add(dedy, degy, ia);
487 atomic_add(dedz, degz, ia);
488 if CONSTEXPR (do_v) {
489 vxx = xr * dedx;
490 vyx = yr * dedx;
491 vzx = zr * dedx;
492 vyy = yr * dedy;
493 vzy = zr * dedy;
494 vzz = zr * dedz;
495 }
496 }
497}
498}
void atomic_add(T value, T *buffer, size_t offset=0)
Definition: acc/adddef.h:14
#define SEQ_CUDA
Definition: acc/seqdef.h:12
#define image(x, y, z)
Definition: image.h:210
#define TINKER_IMAGE_PARAMS
Definition: box.h:110
#define restrict
Definition: macro.h:51
#define CONSTEXPR
Definition: macro.h:61
double * mass
Atomic mass.
real * z
Current coordinates used in energy evaluation and neighbor lists.
real * x
Number of the trajectory frames.
real * y
Current coordinates used in energy evaluation and neighbor lists.
int * ipfix
real * zpfix
grad_prec * degy
real(* tfix)[3]
grad_prec * degx
int(* igfix)[2]
Group numbers defining each group distance restraint.
real(* dfix)[3]
real * ypfix
int(* iafix)[3]
int(* itfix)[4]
int(* kpfix)[3]
real(* pfix)[2]
grad_prec * degz
real(* gfix)[3]
Force constant and target range for each group distance.
int(* idfix)[2]
real(* afix)[3]
real * xpfix
constexpr real _1radian
Definition: const.h:17
constexpr real radian
Definition: const.h:14
float real
Definition: precision.h:80
fixed grad_prec
Definition: precision.h:103
Definition: testrt.h:9
__device__ void dk_geom_group(real &__restrict__ e, real &__restrict__ vxx, real &__restrict__ vyx, real &__restrict__ vzx, real &__restrict__ vyy, real &__restrict__ vzy, real &__restrict__ vzz, grad_prec *__restrict__ degx, grad_prec *__restrict__ degy, grad_prec *__restrict__ degz, int i, const int(*__restrict__ igfix)[2], const real(*__restrict__ gfix)[3], const real *__restrict__ x, const real *__restrict__ y, const real *__restrict__ z, const double *__restrict__ mass, const int *__restrict__ molec, const int(*__restrict__ igrp)[2], const int *__restrict__ kgrp, const double *__restrict__ grpmass, BoxShape box_shape, real3 lvec1, real3 lvec2, real3 lvec3, real3 recipa, real3 recipb, real3 recipc)
Definition: geom.h:12
__device__ void dk_geom_angle(real &__restrict__ e, real &__restrict__ vxx, real &__restrict__ vyx, real &__restrict__ vzx, real &__restrict__ vyy, real &__restrict__ vzy, real &__restrict__ vzz, grad_prec *__restrict__ degx, grad_prec *__restrict__ degy, grad_prec *__restrict__ degz, int i, const int(*__restrict__ iafix)[3], const real(*__restrict__ afix)[3], const real *__restrict__ x, const real *__restrict__ y, const real *__restrict__ z)
Definition: geom.h:192
__device__ void dk_geom_distance(real &__restrict__ e, real &__restrict__ vxx, real &__restrict__ vyx, real &__restrict__ vzx, real &__restrict__ vyy, real &__restrict__ vzy, real &__restrict__ vzz, grad_prec *__restrict__ degx, grad_prec *__restrict__ degy, grad_prec *__restrict__ degz, int i, const int(*__restrict__ idfix)[2], const real(*__restrict__ dfix)[3], const real *__restrict__ x, const real *__restrict__ y, const real *__restrict__ z, const int *__restrict__ molec, BoxShape box_shape, real3 lvec1, real3 lvec2, real3 lvec3, real3 recipa, real3 recipb, real3 recipc)
Definition: geom.h:130
__device__ void dk_geom_position(real &__restrict__ e, real &__restrict__ vxx, real &__restrict__ vyx, real &__restrict__ vzx, real &__restrict__ vyy, real &__restrict__ vzy, real &__restrict__ vzz, grad_prec *__restrict__ degx, grad_prec *__restrict__ degy, grad_prec *__restrict__ degz, int i, const int *__restrict__ ipfix, const int(*__restrict__ kpfix)[3], const real *__restrict__ xpfix, const real *__restrict__ ypfix, const real *__restrict__ zpfix, const real(*__restrict__ pfix)[2], const real *__restrict__ x, const real *__restrict__ y, const real *__restrict__ z, BoxShape box_shape, real3 lvec1, real3 lvec2, real3 lvec3, real3 recipa, real3 recipb, real3 recipc)
Definition: geom.h:447
__device__ void dk_geom_torsion(real &__restrict__ e, real &__restrict__ vxx, real &__restrict__ vyx, real &__restrict__ vzx, real &__restrict__ vyy, real &__restrict__ vzy, real &__restrict__ vzz, grad_prec *__restrict__ degx, grad_prec *__restrict__ degy, grad_prec *__restrict__ degz, int i, const int(*__restrict__ itfix)[4], const real(*__restrict__ tfix)[3], const real *__restrict__ x, const real *__restrict__ y, const real *__restrict__ z)
Definition: geom.h:292