An attempt of adding hard constraints plugin to the simulator. It seems it ...

Samo Penic

2019-02-27 51b4f09bca802c4941700c78c5733ce968ddc7c6

An attempt of adding hard constraints plugin to the simulator. It seems it does not work yet

 src/Makefile.am

@@ -20,3 +20,4 @@
pkginclude_HEADERS=general.h vertex.h bond.h triangle.h cell.h vesicle.h initial_distribution.h io.h dumpstate.h frame.h energy.h timestep.h vertexmove.h bondflip.h poly.h stats.h sh.h shcomplex.h constvol.h snapshot.h restore.h cluster.h plugins.h
#demoplugin_la_SOURCES= plugins/demo_plugin.c
#demoplugin_la_LDFLAGS = -module -avoid-version -export-dynamic
SUBDIRS=plugins

 src/general.h

@@ -306,11 +306,13 @@

/* plugins */
typedef struct {
    void (*at_start)(void);
    void (*after_vesicle_init)(void *vesicle);
    void (*at_start)(int argc, char **argv);
    void *(*after_vesicle_init)(void *vesicle);
    ts_bool (*vm_hard_constraint)(void *vesicle, ts_vertex *vtx, ts_vertex *odl_vtx);
    ts_double (*vm_energy_before)(void *vesicle, ts_vertex *vtx);
    ts_double (*vm_energy_after)(void *vesicle, ts_vertex *vtx);
    ts_double (*vm_energy_before_prepare)(void *vesicle, ts_vertex *vtx);
    ts_double (*vm_energy_before_execute)(void *vesicle, ts_vertex *vtx);
    ts_double (*vm_energy_after_prepare)(void *vesicle, ts_vertex *vtx);
    ts_double (*vm_energy_after_execute)(void *vesicle, ts_vertex *vtx);
    ts_double (*vm_new_state_accepted)(void *vesicle, ts_vertex *vtx);
    ts_double (*vm_new_state_rejected)(void *vesicle, ts_vertex *vtx);
    void (*cleanup)(void);
@@ -365,7 +367,7 @@
    ts_double nucleus_center[3];
    ts_double area;
    ts_confinement_plane confinement_plane;
    ts_plugin_list plist;
    ts_plugin_list *plist;
} ts_vesicle;



 src/main.c

@@ -27,10 +27,13 @@
#include <string.h>
int main(int argv, char *argc[]){
    ts_char *plugin0 = (ts_char *)calloc(255,sizeof(ts_char));
    ts_char *plugin1 = (ts_char *)calloc(255,sizeof(ts_char));
    strcpy(plugin0,"/home/samo/programiranje/trisurf-ng/src/plugins/.libs/demoplugin.so");
    ts_char **plugins=(ts_char **)calloc(1,sizeof(ts_char *));
    strcpy(plugin1,"/home/samo/programiranje/trisurf-ng/src/plugins/.libs/vmdefaulthardconstraints.so");
    ts_char **plugins=(ts_char **)calloc(2,sizeof(ts_char *));
    *plugins=plugin0;
    ts_plugin_list *plist=init_plugin_list(plugins,1);
    plugins[1]=plugin1;
    ts_plugin_list *plist=init_plugin_list(plugins,2);
    ts_fprintf(stdout, "TRISURF in PRVI PLUGIN %s\n", plist->plugin[0]->details->name);
    ts_vesicle *vesicle;
    ts_tape *tape;
@@ -46,15 +49,18 @@
    ts_fprintf(stdout,"Programming done by: Samo Penic and Miha Fosnaric\n");
    ts_fprintf(stdout,"Released under terms of GPLv3\n");
    ts_fprintf(stdout,"Starting program...\n\n");
//    vesicle = parseDump("timestep_000000.vtu");
//        run_simulation(vesicle, vesicle->tape->mcsweeps, vesicle->tape->inititer, vesicle->tape->iterations, 1);
	
/* Entry point for plugin function at_start() */
    int i;
    for(i=0; i<plist->n; i++){
        plist->plugin[i]->function->at_start(argv, argc);
    }

    if(command_line_args.dump_from_vtk[0]!=0){
        ts_fprintf(stdout,"************************************************\n");
        ts_fprintf(stdout,"**** Restoring vesicle from VTK points list ****\n");
        ts_fprintf(stdout,"************************************************\n\n");
        vesicle = parseDump(command_line_args.dump_from_vtk);
//        write_vertex_xml_file(vesicle,9999); // here you can test if restoration and rewritting results in the same dump file. Only the date od creation of dump file must differ.
        tape = vesicle->tape;
        int arguments_no;
        FILE *fd=fopen(".status","r");
@@ -68,9 +74,6 @@
        }
        else
            ts_fprintf(stdout,"No .status file. The iteration count will start from 0\n");
/* Here you should read new tape file, reassign some values in vertex from the tape and assign read tape to vesicle->tape */
//        tape=parsetape(command_line_args.tape_fullfilename);
  //      vesicle=vtk2vesicle(command_line_args.dump_from_vtk,tape);
    }
    else if(command_line_args.force_from_tape){
        ts_fprintf(stdout,"************************************************\n");
@@ -132,11 +135,13 @@
        }
        }
    }
            //printf("nucleus coords: %.17e %.17e %.17e\n",vesicle->nucleus_center[0], vesicle->nucleus_center[1], vesicle->nucleus_center[2]);
//    write_vertex_xml_file(vesicle,0);
//    exit(1);

            //write_vertex_xml_file(vesicle,1000);
    vesicle->plist=plist;
/* Entry point for plugin after_vesicle_init() function */
    for(i=0; i<plist->n; i++){
        vesicle = plist->plugin[i]->function->after_vesicle_init(vesicle);
    }
	
    run_simulation(vesicle, tape->mcsweeps, tape->inititer, tape->iterations, start_iteration);
    write_master_xml_file(command_line_args.output_fullfilename);
    write_dout_fcompat_file(vesicle,"dout");

 src/plugins.c

@@ -110,8 +110,10 @@
    plugin->function->at_start = dlsym(plugin->libhandle, "at_start");
    plugin->function->after_vesicle_init = dlsym(plugin->libhandle, "after_vesicle_init");
    plugin->function->vm_hard_constraint = dlsym(plugin->libhandle, "vm_hard_constraint");
    plugin->function->vm_energy_before = dlsym(plugin->libhandle, "vm_energy_before");
    plugin->function->vm_energy_after = dlsym(plugin->libhandle, "vm_energy_after");
    plugin->function->vm_energy_before_prepare = dlsym(plugin->libhandle, "vm_energy_before_prepare");
    plugin->function->vm_energy_before_execute = dlsym(plugin->libhandle, "vm_energy_before_execute");
    plugin->function->vm_energy_after_prepare = dlsym(plugin->libhandle, "vm_energy_after_prepare");
    plugin->function->vm_energy_after_execute = dlsym(plugin->libhandle, "vm_energy_after_execute");
    plugin->function->vm_new_state_rejected = dlsym(plugin->libhandle, "vm_new_state_rejected");
    plugin->function->vm_new_state_accepted = dlsym(plugin->libhandle, "vm_new_state_accepted");

@@ -125,6 +127,7 @@
    for(int i=0;i<number_of_plugins;i++){
        plist->plugin[i]=init_plugin(plugin_filenames[i]);
    }
    plist->n=number_of_plugins;
    return plist;
}    


 src/plugins/.Makefile.am.swp

Binary files differ

 src/plugins/.libs/demo_plugin.o

Binary files differ

 src/plugins/.libs/demoplugin.a

Binary files differ

 src/plugins/.libs/demoplugin.so

Binary files differ

 src/plugins/Makefile

@@ -137,6 +137,14 @@
demoplugin_la_LINK = $(LIBTOOL) $(AM_V_lt) --tag=CC $(AM_LIBTOOLFLAGS) \
    $(LIBTOOLFLAGS) --mode=link $(CCLD) $(AM_CFLAGS) $(CFLAGS) \
    $(demoplugin_la_LDFLAGS) $(LDFLAGS) -o $@
vmdefaulthardconstraints_la_LIBADD =
am_vmdefaulthardconstraints_la_OBJECTS = default_hard_constraints.lo
vmdefaulthardconstraints_la_OBJECTS =  \
    $(am_vmdefaulthardconstraints_la_OBJECTS)
vmdefaulthardconstraints_la_LINK = $(LIBTOOL) $(AM_V_lt) --tag=CC \
    $(AM_LIBTOOLFLAGS) $(LIBTOOLFLAGS) --mode=link $(CCLD) \
    $(AM_CFLAGS) $(CFLAGS) $(vmdefaulthardconstraints_la_LDFLAGS) \
    $(LDFLAGS) -o $@
AM_V_P = $(am__v_P_$(V))
am__v_P_ = $(am__v_P_$(AM_DEFAULT_VERBOSITY))
am__v_P_0 = false
@@ -171,8 +179,10 @@
am__v_CCLD_ = $(am__v_CCLD_$(AM_DEFAULT_VERBOSITY))
am__v_CCLD_0 = @echo "  CCLD    " $@;
am__v_CCLD_1 = 
SOURCES = $(demoplugin_la_SOURCES)
DIST_SOURCES = $(demoplugin_la_SOURCES)
SOURCES = $(demoplugin_la_SOURCES) \
    $(vmdefaulthardconstraints_la_SOURCES)
DIST_SOURCES = $(demoplugin_la_SOURCES) \
    $(vmdefaulthardconstraints_la_SOURCES)
am__can_run_installinfo = \
  case $$AM_UPDATE_INFO_DIR in \
    n|no|NO) false;; \
@@ -321,9 +331,11 @@
top_builddir = ../..
top_srcdir = ../..
AM_CFLAGS = -Wall -Werror -DTS_VERSION=\"$(GITVERSION)\" -fgnu89-inline -Wno-error=date-time -I..
lib_LTLIBRARIES = demoplugin.la
lib_LTLIBRARIES = demoplugin.la vmdefaulthardconstraints.la
demoplugin_la_SOURCES = demo_plugin.c
demoplugin_la_LDFLAGS = -module -avoid-version -export-dynamic
vmdefaulthardconstraints_la_SOURCES = default_hard_constraints.c
vmdefaulthardconstraints_la_LDFLAGS = -module -avoid-version -export-dynamic
all: all-am

.SUFFIXES:
@@ -396,12 +408,16 @@
demoplugin.la: $(demoplugin_la_OBJECTS) $(demoplugin_la_DEPENDENCIES) $(EXTRA_demoplugin_la_DEPENDENCIES) 
    $(AM_V_CCLD)$(demoplugin_la_LINK) -rpath $(libdir) $(demoplugin_la_OBJECTS) $(demoplugin_la_LIBADD) $(LIBS)

vmdefaulthardconstraints.la: $(vmdefaulthardconstraints_la_OBJECTS) $(vmdefaulthardconstraints_la_DEPENDENCIES) $(EXTRA_vmdefaulthardconstraints_la_DEPENDENCIES) 
    $(AM_V_CCLD)$(vmdefaulthardconstraints_la_LINK) -rpath $(libdir) $(vmdefaulthardconstraints_la_OBJECTS) $(vmdefaulthardconstraints_la_LIBADD) $(LIBS)

mostlyclean-compile:
    -rm -f *.$(OBJEXT)

distclean-compile:
    -rm -f *.tab.c

include ./$(DEPDIR)/default_hard_constraints.Plo
include ./$(DEPDIR)/demo_plugin.Plo

.c.o:

 src/plugins/Makefile.am

@@ -1,5 +1,7 @@
AM_CFLAGS = -Wall -Werror -DTS_VERSION=\"$(GITVERSION)\" -fgnu89-inline -Wno-error=date-time -I..
lib_LTLIBRARIES= demoplugin.la
lib_LTLIBRARIES= demoplugin.la vmdefaulthardconstraints.la
demoplugin_la_SOURCES= demo_plugin.c
demoplugin_la_LDFLAGS = -module -avoid-version -export-dynamic

vmdefaulthardconstraints_la_SOURCES = default_hard_constraints.c
vmdefaulthardconstraints_la_LDFLAGS = -module -avoid-version -export-dynamic

 src/plugins/Makefile.in

@@ -137,6 +137,14 @@
demoplugin_la_LINK = $(LIBTOOL) $(AM_V_lt) --tag=CC $(AM_LIBTOOLFLAGS) \
    $(LIBTOOLFLAGS) --mode=link $(CCLD) $(AM_CFLAGS) $(CFLAGS) \
    $(demoplugin_la_LDFLAGS) $(LDFLAGS) -o $@
vmdefaulthardconstraints_la_LIBADD =
am_vmdefaulthardconstraints_la_OBJECTS = default_hard_constraints.lo
vmdefaulthardconstraints_la_OBJECTS =  \
    $(am_vmdefaulthardconstraints_la_OBJECTS)
vmdefaulthardconstraints_la_LINK = $(LIBTOOL) $(AM_V_lt) --tag=CC \
    $(AM_LIBTOOLFLAGS) $(LIBTOOLFLAGS) --mode=link $(CCLD) \
    $(AM_CFLAGS) $(CFLAGS) $(vmdefaulthardconstraints_la_LDFLAGS) \
    $(LDFLAGS) -o $@
AM_V_P = $(am__v_P_@AM_V@)
am__v_P_ = $(am__v_P_@AM_DEFAULT_V@)
am__v_P_0 = false
@@ -171,8 +179,10 @@
am__v_CCLD_ = $(am__v_CCLD_@AM_DEFAULT_V@)
am__v_CCLD_0 = @echo "  CCLD    " $@;
am__v_CCLD_1 = 
SOURCES = $(demoplugin_la_SOURCES)
DIST_SOURCES = $(demoplugin_la_SOURCES)
SOURCES = $(demoplugin_la_SOURCES) \
    $(vmdefaulthardconstraints_la_SOURCES)
DIST_SOURCES = $(demoplugin_la_SOURCES) \
    $(vmdefaulthardconstraints_la_SOURCES)
am__can_run_installinfo = \
  case $$AM_UPDATE_INFO_DIR in \
    n|no|NO) false;; \
@@ -321,9 +331,11 @@
top_builddir = @top_builddir@
top_srcdir = @top_srcdir@
AM_CFLAGS = -Wall -Werror -DTS_VERSION=\"$(GITVERSION)\" -fgnu89-inline -Wno-error=date-time -I..
lib_LTLIBRARIES = demoplugin.la
lib_LTLIBRARIES = demoplugin.la vmdefaulthardconstraints.la
demoplugin_la_SOURCES = demo_plugin.c
demoplugin_la_LDFLAGS = -module -avoid-version -export-dynamic
vmdefaulthardconstraints_la_SOURCES = default_hard_constraints.c
vmdefaulthardconstraints_la_LDFLAGS = -module -avoid-version -export-dynamic
all: all-am

.SUFFIXES:
@@ -396,12 +408,16 @@
demoplugin.la: $(demoplugin_la_OBJECTS) $(demoplugin_la_DEPENDENCIES) $(EXTRA_demoplugin_la_DEPENDENCIES) 
    $(AM_V_CCLD)$(demoplugin_la_LINK) -rpath $(libdir) $(demoplugin_la_OBJECTS) $(demoplugin_la_LIBADD) $(LIBS)

vmdefaulthardconstraints.la: $(vmdefaulthardconstraints_la_OBJECTS) $(vmdefaulthardconstraints_la_DEPENDENCIES) $(EXTRA_vmdefaulthardconstraints_la_DEPENDENCIES) 
    $(AM_V_CCLD)$(vmdefaulthardconstraints_la_LINK) -rpath $(libdir) $(vmdefaulthardconstraints_la_OBJECTS) $(vmdefaulthardconstraints_la_LIBADD) $(LIBS)

mostlyclean-compile:
    -rm -f *.$(OBJEXT)

distclean-compile:
    -rm -f *.tab.c

@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/default_hard_constraints.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/demo_plugin.Plo@am__quote@

.c.o:

 src/plugins/default_hard_constraints.c

New file
@@ -0,0 +1,60 @@
#include "plugins.h"
#include <stdlib.h>
#include "general.h"
#include "vertex.h"
#include "cell.h"
char plugin_name[] = "Default hard constraints for vertex move";
char plugin_description[]= "Checks proximity with neighbors and with foreign vertices";
char plugin_author[] = "SAMO PENIC";

ts_plugin_details *init (){
//    ts_fprintf(stdout,"Hello. Plugin %s is initiating. This will load the details section of the plugin\n", plugin_name);

    ts_plugin_details *details=(ts_plugin_details *)calloc(1,sizeof(ts_plugin_details));
    details->name = plugin_name;
    return details;	
}

void at_start(int argc, char **argv){
}

void *after_vesicle_init(ts_vesicle *vesicle){
    ts_fprintf(stdout,"Plugin \"%s\" is loaded.\n", plugin_name);
    return vesicle;
}

ts_bool vm_hard_constraint(ts_vesicle *vesicle, ts_vertex *vtx, ts_vertex *ovtx){
    ts_uint i;
    ts_uint cellidx;
    ts_double dist;
    ts_bool retval;

//distance with neighbours check
    for(i=0;i<vtx->neigh_no;i++){
        dist=vtx_distance_sq(vtx,vtx->neigh[i]);
        if(dist<1.0 || dist>vesicle->dmax) {
            return TS_FAIL;
        }
    }
// Distance with grafted poly-vertex check:	
    if(vtx->grafted_poly!=NULL){
        dist=vtx_distance_sq(vtx,vtx->grafted_poly->vlist->vtx[0]);
        if(dist<1.0 || dist>vesicle->dmax) {
        return TS_FAIL;
        }
    }

//self avoidance check with distant vertices
    cellidx=vertex_self_avoidance(vesicle, vtx);
    //check occupation number
    retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,vtx);
        if(retval==TS_FAIL){
           return TS_FAIL;
        } 
   
    return TS_SUCCESS;
}

void cleanup(){
    ts_fprintf(stdout,"Goodbye from plugin %s. This functions clears what would be created in init...\n",plugin_name);
}

 src/plugins/demo_plugin.c

@@ -6,13 +6,29 @@
char plugin_author[] = "SAMO PENIC";

ts_plugin_details *init (){
    ts_fprintf(stdout,"Hello. Plugin %s is initiating.\n\nThis will load the details section of the plugin\n", plugin_name);
//    ts_fprintf(stdout,"Hello. Plugin %s is initiating. This will load the details section of the plugin\n", plugin_name);

    ts_plugin_details *details=(ts_plugin_details *)calloc(1,sizeof(ts_plugin_details));
    details->name = plugin_name;
    return details;    
}


void at_start(int argc, char **argv){
    ts_fprintf(stdout, "This will be run at the start of the program\n");
    ts_fprintf(stdout, "This plugin doesn't do anything useful :/. at_start function will be modified in the future!\n");
}


void *after_vesicle_init(ts_vesicle *vesicle){
    ts_fprintf(stdout, "This function can modify vesicle before the start of the simulation. It has to return pointer to vesicle!!!\n");
    return vesicle;
}

ts_bool vm_hard_constraint(ts_vesicle *vesicle, ts_vertex *vtx, ts_vertex *ovtx){
    return TS_SUCCESS;
}

void cleanup(){
    ts_fprintf(stdout,"Goodbye from plugin %s. This functions clears what would be created in init...\n",plugin_name);
}

 src/plugins/demo_plugin.o

Binary files differ

 src/tape

@@ -12,9 +12,9 @@

# Pressure calculations
# (pswitch=1: calc. p*dV energy contribution)
pswitch = 0
pswitch = 1
# pressure difference: p_inside - p_outside (in units kT/l_min^3):
pressure=10.0
pressure=-100.0

#Constant volume constraint (0 disable constant volume, 1 enable wiht additional vertex move, 2 enable with epsvol)
constvolswitch=0

 src/vertexmove.c

@@ -9,16 +9,15 @@
#include "energy.h"
#include "timestep.h"
#include "cell.h"
//#include "io.h"
#include "io.h"
#include<stdio.h>
#include "vertexmove.h"
#include <string.h>
#include "constvol.h"
#include "plugins.h"

ts_bool single_verticle_timestep(ts_vesicle *vesicle,ts_vertex *vtx,ts_double *rn){
    ts_uint i;
    ts_double dist;
    ts_bool retval; 
    ts_uint cellidx; 
    ts_double delta_energy, delta_energy_cv,oenergy,dvol=0.0, darea=0.0, dstretchenergy=0.0;
@@ -26,21 +25,6 @@
    //This will hold all the information of vtx and its neighbours
    ts_vertex backupvtx[20], *constvol_vtx_moved=NULL, *constvol_vtx_backup=NULL;
    memcpy((void *)&backupvtx[0],(void *)vtx,sizeof(ts_vertex));

    //Some stupid tests for debugging cell occupation!
/*         cellidx=vertex_self_avoidance(vesicle, vtx);
    if(vesicle->clist->cell[cellidx]==vtx->cell){
        fprintf(stderr,"Idx match!\n");
    } else {
        fprintf(stderr,"***** Idx don't match!\n");
        fatal("ENding.",1);
    }
*/

        //temporarly moving the vertex
//    vtx->x=vtx->x+vesicle->stepsize*(2.0*rn[0]-1.0);
//        vtx->y=vtx->y+vesicle->stepsize*(2.0*rn[1]-1.0);
//        vtx->z=vtx->z+vesicle->stepsize*(2.0*rn[2]-1.0);

//random move in a sphere with radius stepsize:
    r=vesicle->stepsize*rn[0];
@@ -51,24 +35,6 @@
    vtx->y=vtx->y+r*sintheta*sin(phi);
    vtx->z=vtx->z+r*costheta;


//distance with neighbours check
    for(i=0;i<vtx->neigh_no;i++){
        dist=vtx_distance_sq(vtx,vtx->neigh[i]);
        if(dist<1.0 || dist>vesicle->dmax) {
            vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex));
            return TS_FAIL;
        }
    }

// Distance with grafted poly-vertex check:	
    if(vtx->grafted_poly!=NULL){
        dist=vtx_distance_sq(vtx,vtx->grafted_poly->vlist->vtx[0]);
        if(dist<1.0 || dist>vesicle->dmax) {
        vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex));
        return TS_FAIL;
        }
    }

// TODO: Maybe faster if checks only nucleus-neighboring cells
// Nucleus penetration check:
@@ -97,16 +63,16 @@
        }
    }
//#undef SQ
//self avoidance check with distant vertices
    cellidx=vertex_self_avoidance(vesicle, vtx);
    //check occupation number
    retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,vtx);

    if(retval==TS_FAIL){
        vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex));
        return TS_FAIL;
    } 
   
/* Entry point for plugin vm_hard_constraint() function */
    for(i=0; i<vesicle->plist->n;i++){
        retval = vesicle->plist->plugin[i]->function->vm_hard_constraint(vesicle,vtx,&backupvtx[0]);
        if(retval==TS_FAIL){
            vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex));
            return TS_FAIL;
        }
    }

 
//if all the tests are successful, then energy for vtx and neighbours is calculated
    for(i=0;i<vtx->neigh_no;i++){
@@ -280,6 +246,7 @@
    //accepted    
 //   fprintf(stderr,"MC accepted\n");
//    oldcellidx=vertex_self_avoidance(vesicle, &backupvtx[0]);
    cellidx=vertex_self_avoidance(vesicle, vtx);
    if(vtx->cell!=vesicle->clist->cell[cellidx]){
        retval=cell_add_vertex(vesicle->clist->cell[cellidx],vtx);
//        if(retval==TS_SUCCESS) cell_remove_vertex(vesicle->clist->cell[oldcellidx],vtx);