From 19c1fd84d1fd75e9bfd55feb670fda87a2e02086 Mon Sep 17 00:00:00 2001
From: Samo Penic <samo.penic@gmail.com>
Date: Fri, 08 Feb 2019 08:32:29 +0000
Subject: [PATCH] Merge branch 'master' of ssh://git.penic.eu:29418/trisurf-ng

---
 src/constvol.c |   47 ++++++++++++++++++++++++++++++-----------------
 1 files changed, 30 insertions(+), 17 deletions(-)

diff --git a/src/constvol.c b/src/constvol.c
index be69edf..63c07e0 100644
--- a/src/constvol.c
+++ b/src/constvol.c
@@ -1,3 +1,4 @@
+/* vim: set ts=4 sts=4 sw=4 noet : */
 #include<stdlib.h>
 #include<stdio.h>
 #include<string.h>
@@ -12,7 +13,7 @@
 ts_bool constvolume(ts_vesicle *vesicle, ts_vertex *vtx_avoid, ts_double Vol, ts_double *retEnergy, ts_vertex **vtx_moved_retval, ts_vertex **vtx_backup){
     ts_vertex *vtx_moved;
     ts_uint vtxind,i,j;
-    ts_uint Ntries=3;
+    ts_uint Ntries=100;
 	ts_vertex *backupvtx;
     ts_double Rv, dh, dvol, volFirst, voldiff, oenergy,delta_energy;
     backupvtx=(ts_vertex *)calloc(sizeof(ts_vertex),10);
@@ -20,13 +21,22 @@
     for(i=0;i<Ntries;i++){
         vtxind=rand() % vesicle->vlist->n;
         vtx_moved=vesicle->vlist->vtx[vtxind];
-        /* chosen vertex must not be a nearest neighbour. TODO: probably must
+        /* chosen vertex must not be a nearest neighbour. WASTODO: probably must.
+         * DONE: solved below, when using different algorithm.
          * extend search in case of bondflip */
-        if(vtx_moved==vtx_avoid) continue;
+/*        if(vtx_moved==vtx_avoid) continue;
         for(j=0;j<vtx_moved->neigh_no;j++){
             if(vtx_moved->neigh[j]==vtx_avoid) continue;
         }
-         
+*/
+/* different check of nearest neighbour (and second nearest neighbour) check.
+ * Checking the distance between vertex and vertex to be moved to assure
+ * constant volume. Solves upper todo problem. */
+
+        if(vtx_distance_sq(vtx_moved,vtx_avoid)<16.0*vesicle->dmax){
+            continue;
+        }
+        
 	    memcpy((void *)&backupvtx[0],(void *)vtx_moved,sizeof(ts_vertex));
 
         //move vertex in specified direction. first try, test move!
@@ -36,16 +46,16 @@
 	    vtx_moved->y=vtx_moved->y*(1.0-dh/Rv);
 	    vtx_moved->z=vtx_moved->z*(1.0-dh/Rv);
 
+
+//SKIPPING FIRST CHECK of CONSTRAINTS. This is not a final move.
         //check for constraints
-          if(constvolConstraintCheck(vesicle, vtx_moved)==TS_FAIL){
+/*          if(constvolConstraintCheck(vesicle, vtx_moved)==TS_FAIL){
 		    vtx_moved=memcpy((void *)vtx_moved,(void *)&backupvtx[0],sizeof(ts_vertex));
-            continue;
-        }
+//            continue;
+                break;
+        }  */
         // All checks OK!
 
-        for(j=0;j<vtx_moved->neigh_no;j++){
-        	memcpy((void *)&backupvtx[j+1],(void *)vtx_moved->neigh[j],sizeof(ts_vertex));
-	    }
         dvol=0.0;
         for(j=0;j<vtx_moved->tristar_no;j++){
             dvol-=vtx_moved->tristar[j]->volume;
@@ -86,13 +96,11 @@
 	    vtx_moved->z=vtx_moved->z*(1-dh/Rv);
         //check for constraints
         if(constvolConstraintCheck(vesicle, vtx_moved)==TS_FAIL){
-	        for(j=0;j<vtx_moved->neigh_no;j++){
-        	    memcpy((void *)vtx_moved->neigh[j],(void *)&backupvtx[j+1],sizeof(ts_vertex));
-	        }
 	        vtx_moved=memcpy((void *)vtx_moved,(void *)&backupvtx[0],sizeof(ts_vertex));
             //also, restore normals
             for(j=0;j<vtx_moved->tristar_no;j++) triangle_normal_vector(vtx_moved->tristar[j]);
-            continue;
+//            continue;
+                break;
         }
 
         dvol=volFirst;
@@ -104,6 +112,10 @@
         if(fabs(voldiff)/vesicle->volume < vesicle->tape->constvolprecision){
             //calculate energy, return change in energy...
 //            fprintf(stderr, "Constvol success! %e\n",voldiff);
+            for(j=0;j<vtx_moved->neigh_no;j++){
+        	memcpy((void *)&backupvtx[j+1],(void *)vtx_moved->neigh[j],sizeof(ts_vertex));
+    	    }
+
             oenergy=vtx_moved->energy;
             energy_vertex(vtx_moved);
             delta_energy=vtx_moved->xk*(vtx_moved->energy - oenergy);
@@ -159,10 +171,11 @@
 ts_bool constvolumerestore(ts_vertex *vtx_moved,ts_vertex *vtx_backup){
     ts_uint j;
 	 memcpy((void *)vtx_moved,(void *)&vtx_backup[0],sizeof(ts_vertex));
-     for(j=0;j<vtx_moved->neigh_no;j++){
-        	    memcpy((void *)vtx_moved->neigh[j],(void *)&vtx_backup[j+1],sizeof(ts_vertex));
-	}
     for(j=0;j<vtx_moved->tristar_no;j++) triangle_normal_vector(vtx_moved->tristar[j]);
+     for(j=0;j<vtx_moved->neigh_no;j++){
+        	   // memcpy((void *)vtx_moved->neigh[j],(void *)&vtx_backup[j+1],sizeof(ts_vertex));
+                 energy_vertex(vtx_moved->neigh[j]);
+	}
 
     free(vtx_backup);
     return TS_SUCCESS;

--
Gitblit v1.9.3