nyereg pont + poincare-hopf

epgpu.cu

@@ -148,10 +148,38 @@ __device__ inline vec3 lineintersect(const vec3& S, const vec3& A, const vec3& B
     return A + (dot(S - A, B - A) * (B - A)) / dot(B - A, B - A); 
     return A + (dot(S - A, B - A) * (B - A)) / dot(B - A, B - A); 
 }
 }
 
 
+__device__ inline bool nyereg_ell(const vec3& S, const vec3& X1, const vec3& X2, const vec3& A, const vec3& B) {
+    vec3 intersect = lineintersect(S, X1, X2);
+    double KAC = dot(X2 - X1, S - X1);
+    double KAB = dot(X2 - X1, X2 - X1);
+    if (0.0 < KAC && KAC < KAB) {
+        vec3 planeN = normalize(intersect - S);
+        double signDistanceA = signeddistance(planeN, intersect, A);
+        double signDistanceB = signeddistance(planeN, intersect, B);
+        if (signDistanceA <= 0.0 && signDistanceB <= 0.0)
+            return true;
+        if (signDistanceA > 0.0 && signDistanceB > 0.0)
+            return true;
+    }
+    return false;
+}
+
 __device__ inline int nyereg_ep(const vec3& S, const vec3& C, const vec3& D) {
 __device__ inline int nyereg_ep(const vec3& S, const vec3& C, const vec3& D) {
     int cnt = 0;
     int cnt = 0;
     vec3 A(0, 0, 0), B(0, 0, 1.0);
     vec3 A(0, 0, 0), B(0, 0, 1.0);
 
 
+    if (nyereg_ell(S, A, B, C, D))
+        cnt++;
+    if (nyereg_ell(S, B, C, D, A))
+        cnt++;
+    if (nyereg_ell(S, C, D, A, B))
+        cnt++;
+    if (nyereg_ell(S, A, D, C, A))
+        cnt++;
+    if (nyereg_ell(S, A, C, B, D))
+        cnt++;
+    if (nyereg_ell(S, D, B, C, A))
+        cnt++;
 
 
     return cnt;
     return cnt;
 }
 }
@@ -173,9 +201,9 @@ __device__ void ABC_oldal(int v, int w, const vec3& C, const vec3& D, char* egys
                 vec3 Sv = K + L*k;
                 vec3 Sv = K + L*k;
                 int S = stabil_ep(Sv, C, D);
                 int S = stabil_ep(Sv, C, D);
                 int U = instabil_ep(Sv, C, D);
                 int U = instabil_ep(Sv, C, D);
-                int H = 0;
+                int H = nyereg_ep(Sv, C, D);;
 
 
-                if (S > 0 && U > 0)
+                if (S > 0 && U > 0 && S + U - H == 2)
                     egysulyi_mtx[pos*16+(S-1)*4+(U-1)] = 1;
                     egysulyi_mtx[pos*16+(S-1)*4+(U-1)] = 1;
             }
             }
         }
         }

tetrarun.py

@@ -40,7 +40,7 @@ def main(argv):
 
 
    res = start_kernel(Cx, Cy, Dx, Dy, Dz, v, w)
    res = start_kernel(Cx, Cy, Dx, Dy, Dz, v, w)
 
 
-   print(res)
+   #print(res)
 
 
 if __name__ == "__main__":
 if __name__ == "__main__":
    main(sys.argv[1:])
    main(sys.argv[1:])
\ No newline at end of file