svn merge -r 12684:12691 https://svn.blender.org/svnroot/bf-blender/trunk/blender
authorDaniel Genrich <daniel.genrich@gmx.net>
Tue, 27 Nov 2007 13:25:07 +0000 (13:25 +0000)
committerDaniel Genrich <daniel.genrich@gmx.net>
Tue, 27 Nov 2007 13:25:07 +0000 (13:25 +0000)
1  2 
source/blender/blenkernel/BKE_cloth.h
source/blender/blenkernel/intern/cloth.c
source/blender/blenkernel/intern/implicit.c
source/blender/blenkernel/intern/modifier.c
source/blender/blenkernel/intern/softbody.c
source/blender/src/buttons_editing.c
source/blender/src/buttons_object.c

index 5a986119b71b4f70b6633293e1927e81043e003b,0000000000000000000000000000000000000000..2bb3e32d6b688861dbed4c357b958ffc60bc294e
mode 100644,000000..100644
--- /dev/null
@@@ -1,275 -1,0 +1,276 @@@
 +/**
 + * BKE_cloth.h
 + *
 + * $Id: BKE_cloth.h,v 1.1 2007/08/01 02:07:27 daniel Exp $
 + *
 + * ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
 + *
 + * This program is free software; you can redistribute it and/or
 + * modify it under the terms of the GNU General Public License
 + * as published by the Free Software Foundation; either version 2
 + * of the License, or (at your option) any later version. The Blender
 + * Foundation also sells licenses for use in proprietary software under
 + * the Blender License.  See http://www.blender.org/BL/ for information
 + * about this.
 + *
 + * This program is distributed in the hope that it will be useful,
 + * but WITHOUT ANY WARRANTY; without even the implied warranty of
 + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 + * GNU General Public License for more details.
 + *
 + * You should have received a copy of the GNU General Public License
 + * along with this program; if not, write to the Free Software Foundation,
 + * Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
 + *
 + * The Original Code is Copyright (C) Blender Foundation.
 + * All rights reserved.
 + *
 + * The Original Code is: all of this file.
 + *
 + * Contributor(s): none yet.
 + *
 + * ***** END GPL/BL DUAL LICENSE BLOCK *****
 + */
 +#ifndef BKE_CLOTH_H
 +#define BKE_CLOTH_H
 +
 +#include "BKE_customdata.h"
 +#include "BLI_linklist.h"
 +#include "BKE_DerivedMesh.h"
 +#include "BKE_object.h"
 +
 +#include "DNA_cloth_types.h"
 +#include "DNA_customdata_types.h"
 +#include "DNA_meshdata_types.h"
 +
 +struct Object;
 +struct Cloth;
 +struct MFace;
 +struct DerivedMesh;
 +struct ClothModifierData;
 +
 +// this is needed for inlining behaviour
 +
 +#ifndef _WIN32
 +#define LINUX
 +#define DO_INLINE inline
 +#else
 +#define DO_INLINE
 +#endif
 +
 +#define CLOTH_MAX_THREAD 2
 +
 +
 +typedef struct ClothVertex {
 +      int     flags;          /* General flags per vertex.            */
 +      float   mass;           /* mass / weight of the vertex          */
 +      float   goal;           /* goal, from SB                        */
 +      float   impulse[3];     /* used in collision.c */
 +      unsigned int impulse_count; /* same as above */
 +      float collball;
++      char octantflag;
 +} ClothVertex;
 +
 +typedef struct ClothSpring {
 +      int     ij;             /* Pij from the paper, one end of the spring.   */
 +      int     kl;             /* Pkl from the paper, one end of the spring.   */
 +      float   restlen;        /* The original length of the spring.   */
 +      int     matrix_index;   /* needed for implicit solver (fast lookup) */
 +      int     type;           /* types defined in BKE_cloth.h ("springType") */
 +      int     flags;          /* defined in BKE_cloth.h, e.g. deactivated due to tearing */
 +      float dfdx[3][3];
 +      float dfdv[3][3];
 +      float f[3];
 +} ClothSpring;
 +
 +typedef struct Cloth {
 +      struct ClothVertex      *verts;                 /* The vertices that represent this cloth. */
 +      struct LinkNode         *springs;               /* The springs connecting the mesh. */
 +      struct BVH              *tree;          /* collision tree for this cloth object */
 +      struct BVH              *selftree;              /* self collision tree for this cloth object */
 +      struct MFace            *mfaces;
 +      struct Implicit_Data    *implicit;      /* our implicit solver connects to this pointer */
 +      struct EdgeHash         *edgehash; /* used for fast checking adjacent points */
 +      unsigned int            numverts;               /* The number of verts == m * n. */
 +      unsigned int            numsprings;             /* The count of springs. */
 +      unsigned int            numfaces;
 +      unsigned int            numothersprings;
 +      unsigned int            numspringssave;
 +      unsigned int            old_solver_type;
 +      float                   (*x)[3]; /* The current position of all vertices.*/
 +      float                   (*xold)[3]; /* The previous position of all vertices.*/
 +      float                   (*current_x)[3]; /* The TEMPORARY current position of all vertices.*/
 +      float                   (*current_xold)[3]; /* The TEMPORARY previous position of all vertices.*/
 +      float                   (*v)[4]; /* the current velocity of all vertices */
 +      float                   (*current_v)[3];
 +      float                   (*xconst)[3];
 +} Cloth;
 +
 +/* goal defines */
 +#define SOFTGOALSNAP  0.999f
 +
 +/* This is approximately the smallest number that can be
 +* represented by a float, given its precision. */
 +#define ALMOST_ZERO           0.000001
 +
 +// some macro enhancements for vector treatment
 +#define VECADDADD(v1,v2,v3)   {*(v1)+= *(v2) + *(v3); *(v1+1)+= *(v2+1) + *(v3+1); *(v1+2)+= *(v2+2) + *(v3+2);}
 +#define VECSUBADD(v1,v2,v3)   {*(v1)-= *(v2) + *(v3); *(v1+1)-= *(v2+1) + *(v3+1); *(v1+2)-= *(v2+2) + *(v3+2);}
 +#define VECADDSUB(v1,v2,v3)   {*(v1)+= *(v2) - *(v3); *(v1+1)+= *(v2+1) - *(v3+1); *(v1+2)+= *(v2+2) - *(v3+2);}
 +#define VECSUBADDSS(v1,v2,aS,v3,bS)   {*(v1)-= *(v2)*aS + *(v3)*bS; *(v1+1)-= *(v2+1)*aS + *(v3+1)*bS; *(v1+2)-= *(v2+2)*aS + *(v3+2)*bS;}
 +#define VECADDSUBSS(v1,v2,aS,v3,bS)   {*(v1)+= *(v2)*aS - *(v3)*bS; *(v1+1)+= *(v2+1)*aS - *(v3+1)*bS; *(v1+2)+= *(v2+2)*aS - *(v3+2)*bS;}
 +#define VECADDSS(v1,v2,aS,v3,bS)      {*(v1)= *(v2)*aS + *(v3)*bS; *(v1+1)= *(v2+1)*aS + *(v3+1)*bS; *(v1+2)= *(v2+2)*aS + *(v3+2)*bS;}
 +#define VECADDS(v1,v2,v3,bS)  {*(v1)= *(v2) + *(v3)*bS; *(v1+1)= *(v2+1) + *(v3+1)*bS; *(v1+2)= *(v2+2) + *(v3+2)*bS;}
 +#define VECSUBMUL(v1,v2,aS)   {*(v1)-= *(v2) * aS; *(v1+1)-= *(v2+1) * aS; *(v1+2)-= *(v2+2) * aS;}
 +#define VECSUBS(v1,v2,v3,bS)  {*(v1)= *(v2) - *(v3)*bS; *(v1+1)= *(v2+1) - *(v3+1)*bS; *(v1+2)= *(v2+2) - *(v3+2)*bS;}
 +#define VECSUBSB(v1,v2, v3,bS)        {*(v1)= (*(v2)- *(v3))*bS; *(v1+1)= (*(v2+1) - *(v3+1))*bS; *(v1+2)= (*(v2+2) - *(v3+2))*bS;}
 +#define VECMULS(v1,aS)        {*(v1)*= aS; *(v1+1)*= aS; *(v1+2)*= *aS;}
 +#define VECADDMUL(v1,v2,aS)   {*(v1)+= *(v2) * aS; *(v1+1)+= *(v2+1) * aS; *(v1+2)+= *(v2+2) * aS;}
 +
 +/* SIMULATION FLAGS: goal flags,.. */
 +/* These are the bits used in SimSettings.flags. */
 +typedef enum
 +{
 +    CLOTH_SIMSETTINGS_FLAG_RESET = ( 1 << 1 ),                // The CM object requires a reinitializaiton.
 +    CLOTH_SIMSETTINGS_FLAG_COLLOBJ = ( 1 << 2 ),      // object is only collision object, no cloth simulation is done
 +    CLOTH_SIMSETTINGS_FLAG_GOAL = ( 1 << 3 ),                 // we have goals enabled
 +    CLOTH_SIMSETTINGS_FLAG_TEARING = ( 1 << 4 ), // true if tearing is enabled
 +    CLOTH_SIMSETTINGS_FLAG_CCACHE_PROTECT = ( 1 << 5 ),
 +    CLOTH_SIMSETTINGS_FLAG_BIG_FORCE = ( 1 << 6 ), // true if we have big spring force for bending
 +    CLOTH_SIMSETTINGS_FLAG_SLEEP = ( 1 << 7 ), // true if we let the cloth go to sleep
 +} CLOTH_SIMSETTINGS_FLAGS;
 +
 +/* SPRING FLAGS */
 +typedef enum
 +{
 +    CLOTH_COLLISIONSETTINGS_FLAG_ENABLED = ( 1 << 1 ),
 +} CLOTH_COLLISIONSETTINGS_FLAGS;
 +
 +/* Spring types as defined in the paper.*/
 +typedef enum
 +{
 +    CLOTH_SPRING_TYPE_STRUCTURAL = 0,
 +    CLOTH_SPRING_TYPE_SHEAR,
 +    CLOTH_SPRING_TYPE_BENDING,
 +} CLOTH_SPRING_TYPES;
 +
 +/* SPRING FLAGS */
 +typedef enum
 +{
 +    CLOTH_SPRING_FLAG_DEACTIVATE = ( 1 << 1 ),
 +    CLOTH_SPRING_FLAG_NEEDED = ( 1 << 2 ), // springs has values to be applied
 +} CLOTH_SPRINGS_FLAGS;
 +
 +/* Bits to or into the ClothVertex.flags. */
 +#define CVERT_FLAG_PINNED     1
 +#define CVERT_FLAG_COLLISION  2
 +
 +
 +// needed for buttons_object.c
 +void cloth_clear_cache(struct Object *ob, struct ClothModifierData *clmd, float framenr);
 +void cloth_free_modifier ( struct ClothModifierData *clmd );
 +
 +// needed for cloth.c
 +void implicit_set_positions ( struct ClothModifierData *clmd );
 +
 +// from cloth.c, needed for modifier.c
 +DerivedMesh *clothModifier_do(struct ClothModifierData *clmd, struct Object *ob, struct DerivedMesh *dm, int useRenderParams, int isFinalCalc);
 +
 +// needed in implicit.c
 +int cloth_bvh_objcollision(struct ClothModifierData *clmd, float step, float prevstep, float dt);
 +
 +////////////////////////////////////////////////
 +
 +
 +/////////////////////////////////////////////////
 +// cloth.c
 +////////////////////////////////////////////////
 +void cloth_free_modifier ( struct ClothModifierData *clmd );
 +void cloth_init ( struct ClothModifierData *clmd );
 +////////////////////////////////////////////////
 +
 +
 +/* Typedefs for function pointers we need for solvers and collision detection. */
 +typedef void ( *CM_COLLISION_SELF ) ( struct ClothModifierData *clmd, int step );
 +// typedef void ( *CM_COLLISION_OBJ ) ( ClothModifierData *clmd, int step, CM_COLLISION_RESPONSE collision_response );
 +
 +
 +/* This enum provides the IDs for our solvers. */
 +// only one available in the moment
 +typedef enum {
 +    CM_IMPLICIT = 0,
 +    CM_VERLET = 1,
 +} CM_SOLVER_ID;
 +
 +
 +/* This structure defines how to call the solver.
 +*/
 +typedef struct
 +{
 +      char            *name;
 +      CM_SOLVER_ID    id;
 +      int     ( *init ) ( struct Object *ob, struct ClothModifierData *clmd );
 +      int     ( *solver ) ( struct Object *ob, float framenr, struct ClothModifierData *clmd, struct ListBase *effectors );
 +      int     ( *free ) ( struct ClothModifierData *clmd );
 +}
 +CM_SOLVER_DEF;
 +
 +
 +/* new C implicit simulator */
 +int implicit_init ( struct Object *ob, struct ClothModifierData *clmd );
 +int implicit_free ( struct ClothModifierData *clmd );
 +int implicit_solver ( struct Object *ob, float frame, struct ClothModifierData *clmd, struct ListBase *effectors );
 +
 +/* explicit verlet simulator */
 +int verlet_init ( struct Object *ob, struct ClothModifierData *clmd );
 +int verlet_free ( struct ClothModifierData *clmd );
 +int verlet_solver ( struct Object *ob, float frame, struct ClothModifierData *clmd, struct ListBase *effectors );
 +
 +
 +/* used for collisions in collision.c */
 +/*
 +typedef struct CollPair
 +{
 +      unsigned int face1; // cloth face
 +      unsigned int face2; // object face
 +      double distance; // magnitude of vector
 +      float normal[3];
 +      float vector[3]; // unnormalized collision vector: p2-p1
 +      float pa[3], pb[3]; // collision point p1 on face1, p2 on face2
 +      int lastsign; // indicates if the distance sign has changed, unused itm
 +      float time; // collision time, from 0 up to 1
 +      unsigned int ap1, ap2, ap3, bp1, bp2, bp3, bp4;
 +      unsigned int pointsb[4];
 +}
 +CollPair;
 +*/
 +
 +/* used for collisions in collision.c */
 +typedef struct EdgeCollPair
 +{
 +      unsigned int p11, p12, p21, p22;
 +      float normal[3];
 +      float vector[3];
 +      float time;
 +      int lastsign;
 +      float pa[3], pb[3]; // collision point p1 on face1, p2 on face2
 +}
 +EdgeCollPair;
 +
 +/* used for collisions in collision.c */
 +typedef struct FaceCollPair
 +{
 +      unsigned int p11, p12, p13, p21;
 +      float normal[3];
 +      float vector[3];
 +      float time;
 +      int lastsign;
 +      float pa[3], pb[3]; // collision point p1 on face1, p2 on face2
 +}
 +FaceCollPair;
 +
 +// function definitions from implicit.c
 +DO_INLINE void mul_fvector_S(float to[3], float from[3], float scalar);
 +
 +#endif
 +
index b997049a637ec47d19e4f1efa8948997bf47da66,0000000000000000000000000000000000000000..6eacb26c315dd4c5c0049b08cc7c31bdf25462ba
mode 100644,000000..100644
--- /dev/null
@@@ -1,1213 -1,0 +1,1224 @@@
- int cloth_build_springs ( Cloth *cloth, DerivedMesh *dm );
 +/*  cloth.c
 +*
 +*
 +* ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
 +*
 +* This program is free software; you can redistribute it and/or
 +* modify it under the terms of the GNU General Public License
 +* as published by the Free Software Foundation; either version 2
 +* of the License, or (at your option) any later version. The Blender
 +* Foundation also sells licenses for use in proprietary software under
 +* the Blender License.  See http://www.blender.org/BL/ for information
 +* about this.
 +*
 +* This program is distributed in the hope that it will be useful,
 +* but WITHOUT ANY WARRANTY; without even the implied warranty of
 +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 +* GNU General Public License for more details.
 +*
 +* You should have received a copy of the GNU General Public License
 +* along with this program; if not, write to the Free Software Foundation,
 +* Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
 +*
 +* The Original Code is Copyright (C) Blender Foundation
 +* All rights reserved.
 +*
 +* The Original Code is: all of this file.
 +*
 +* Contributor(s): none yet.
 +*
 +* ***** END GPL/BL DUAL LICENSE BLOCK *****
 +*/
 +
 +
 +#include <math.h>
 +#include <stdlib.h>
 +#include <string.h>
 +
 +#include "MEM_guardedalloc.h"
 +
 +/* types */
 +#include "DNA_curve_types.h"
 +#include "DNA_cloth_types.h"
 +#include "DNA_object_types.h"
 +#include "DNA_object_force.h"
 +#include "DNA_key_types.h"
 +#include "DNA_lattice_types.h"
 +#include "DNA_mesh_types.h"
 +#include "DNA_meshdata_types.h"
 +#include "DNA_modifier_types.h"
 +#include "DNA_scene_types.h"
 +
 +#include "BLI_blenlib.h"
 +#include "BLI_arithb.h"
 +#include "BLI_edgehash.h"
 +#include "BLI_linklist.h"
 +
 +#include "BKE_curve.h"
 +#include "BKE_cloth.h"
 +#include "BKE_collisions.h"
 +#include "BKE_deform.h"
 +#include "BKE_DerivedMesh.h"
 +#include "BKE_cdderivedmesh.h"
 +#include "BKE_displist.h"
 +#include "BKE_effect.h"
 +#include "BKE_global.h"
 +#include "BKE_key.h"
 +#include "BKE_mesh.h"
 +#include "BKE_modifier.h"
 +#include "BKE_object.h"
 +#include "BKE_pointcache.h"
 +#include "BKE_utildefines.h"
 +
 +#include "BIF_editdeform.h"
 +#include "BIF_editkey.h"
 +#include "DNA_screen_types.h"
 +#include "BSE_headerbuttons.h"
 +#include "BIF_screen.h"
 +#include "BIF_space.h"
 +#include "mydevice.h"
 +
 +#ifdef _WIN32
 +void tstart ( void )
 +{}
 +void tend ( void )
 +{
 +}
 +double tval()
 +{
 +      return 0;
 +}
 +#else
 +#include <sys/time.h>
 +                       static struct timeval _tstart, _tend;
 +       static struct timezone tz;
 +       void tstart ( void )
 +{
 +      gettimeofday ( &_tstart, &tz );
 +}
 +void tend ( void )
 +{
 +      gettimeofday ( &_tend,&tz );
 +}
 +double tval()
 +{
 +      double t1, t2;
 +      t1 = ( double ) _tstart.tv_sec + ( double ) _tstart.tv_usec/ ( 1000*1000 );
 +      t2 = ( double ) _tend.tv_sec + ( double ) _tend.tv_usec/ ( 1000*1000 );
 +      return t2-t1;
 +}
 +#endif
 +
 +/* Our available solvers. */
 +// 255 is the magic reserved number, so NEVER try to put 255 solvers in here!
 +// 254 = MAX!
 +static CM_SOLVER_DEF solvers [] =
 +{
 +      { "Implicit", CM_IMPLICIT, implicit_init, implicit_solver, implicit_free },
 +      // { "Verlet", CM_VERLET, verlet_init, verlet_solver, verlet_free },
 +};
 +
 +/* ********** cloth engine ******* */
 +/* Prototypes for internal functions.
 +*/
 +static void cloth_to_object (Object *ob,  DerivedMesh *dm, ClothModifierData *clmd);
 +static void cloth_from_mesh (Object *ob, ClothModifierData *clmd, DerivedMesh *dm, float framenr);
 +static int cloth_from_object(Object *ob, ClothModifierData *clmd, DerivedMesh *dm, DerivedMesh *olddm, float framenr);
 +
-       clmd->coll_parms->selfepsilon = 0.1;
++int cloth_build_springs ( ClothModifierData *clmd, DerivedMesh *dm );
 +static void cloth_apply_vgroup(ClothModifierData *clmd, DerivedMesh *dm, short vgroup);
 +/******************************************************************************
 +*
 +* External interface called by modifier.c clothModifier functions.
 +*
 +******************************************************************************/
 +/**
 + * cloth_init -  creates a new cloth simulation.
 + *
 + * 1. create object
 + * 2. fill object with standard values or with the GUI settings if given 
 + */
 +void cloth_init (ClothModifierData *clmd)
 +{
 +      /* Initialize our new data structure to reasonable values. */
 +      clmd->sim_parms->gravity [0] = 0.0;
 +      clmd->sim_parms->gravity [1] = 0.0;
 +      clmd->sim_parms->gravity [2] = -9.81;
 +      clmd->sim_parms->structural = 100.0;
 +      clmd->sim_parms->shear = 100.0;
 +      clmd->sim_parms->bending = 1.0;
 +      clmd->sim_parms->Cdis = 5.0;
 +      clmd->sim_parms->Cvi = 1.0;
 +      clmd->sim_parms->mass = 1.0f;
 +      clmd->sim_parms->stepsPerFrame = 5;
 +      clmd->sim_parms->sim_time = 1.0;
 +      clmd->sim_parms->flags = CLOTH_SIMSETTINGS_FLAG_RESET;
 +      clmd->sim_parms->solver_type = 0; 
 +      clmd->sim_parms->preroll = 0;
 +      clmd->sim_parms->maxspringlen = 10;
 +      clmd->coll_parms->self_friction = 5.0;
 +      clmd->coll_parms->friction = 10.0;
 +      clmd->coll_parms->loop_count = 1;
 +      clmd->coll_parms->epsilon = 0.01;
-                               if (!cloth_build_springs (clmd->clothObject, dm) )
++      clmd->coll_parms->selfepsilon = 0.49;
 +      
 +      /* These defaults are copied from softbody.c's
 +      * softbody_calc_forces() function.
 +      */
 +      clmd->sim_parms->eff_force_scale = 1000.0;
 +      clmd->sim_parms->eff_wind_scale = 250.0;
 +
 +      // also from softbodies
 +      clmd->sim_parms->maxgoal = 1.0;
 +      clmd->sim_parms->mingoal = 0.0;
 +      clmd->sim_parms->defgoal = 0.0;
 +      clmd->sim_parms->goalspring = 100.0;
 +      clmd->sim_parms->goalfrict = 0.0;
 +}
 +
 +// unused in the moment, cloth needs quads from mesh
 +DerivedMesh *CDDM_convert_to_triangle(DerivedMesh *dm)
 +{
 +      /*
 +      DerivedMesh *result = NULL;
 +      int i;
 +      int numverts = dm->getNumVerts(dm);
 +      int numedges = dm->getNumEdges(dm);
 +      int numfaces = dm->getNumFaces(dm);
 +
 +      MVert *mvert = CDDM_get_verts(dm);
 +      MEdge *medge = CDDM_get_edges(dm);
 +      MFace *mface = CDDM_get_faces(dm);
 +
 +      MVert *mvert2;
 +      MFace *mface2;
 +      unsigned int numtris=0;
 +      unsigned int numquads=0;
 +      int a = 0;
 +      int random = 0;
 +      int firsttime = 0;
 +      float vec1[3], vec2[3], vec3[3], vec4[3], vec5[3];
 +      float mag1=0, mag2=0;
 +
 +      for(i = 0; i < numfaces; i++)
 +      {
 +      if(mface[i].v4)
 +      numquads++;
 +      else
 +      numtris++;      
 +}
 +
 +      result = CDDM_from_template(dm, numverts, 0, numtris + 2*numquads);
 +
 +      if(!result)
 +      return NULL;
 +
 +      // do verts
 +      mvert2 = CDDM_get_verts(result);
 +      for(a=0; a<numverts; a++) 
 +      {
 +      MVert *inMV;
 +      MVert *mv = &mvert2[a];
 +
 +      inMV = &mvert[a];
 +
 +      DM_copy_vert_data(dm, result, a, a, 1);
 +      *mv = *inMV;
 +}
 +
 +
 +      // do faces
 +      mface2 = CDDM_get_faces(result);
 +      for(a=0, i=0; a<numfaces; a++) 
 +      {
 +      MFace *mf = &mface2[i];
 +      MFace *inMF;
 +      inMF = &mface[a];
 +
 +              
 +              // DM_copy_face_data(dm, result, a, i, 1);
 +
 +              // *mf = *inMF;
 +              
 +
 +      if(mface[a].v4 && random==1)
 +      {
 +      mf->v1 = mface[a].v2;
 +      mf->v2 = mface[a].v3;
 +      mf->v3 = mface[a].v4;
 +}
 +      else
 +      {
 +      mf->v1 = mface[a].v1;
 +      mf->v2 = mface[a].v2;
 +      mf->v3 = mface[a].v3;
 +}
 +
 +      mf->v4 = 0;
 +      mf->flag |= ME_SMOOTH;
 +
 +      test_index_face(mf, NULL, 0, 3);
 +
 +      if(mface[a].v4)
 +      {
 +      MFace *mf2;
 +
 +      i++;
 +
 +      mf2 = &mface2[i];
 +                      
 +                      // DM_copy_face_data(dm, result, a, i, 1);
 +
 +                      // *mf2 = *inMF;
 +                      
 +
 +      if(random==1)
 +      {
 +      mf2->v1 = mface[a].v1;
 +      mf2->v2 = mface[a].v2;
 +      mf2->v3 = mface[a].v4;
 +}
 +      else
 +      {
 +      mf2->v1 = mface[a].v4;
 +      mf2->v2 = mface[a].v1;
 +      mf2->v3 = mface[a].v3;
 +}
 +      mf2->v4 = 0;
 +      mf2->flag |= ME_SMOOTH;
 +
 +      test_index_face(mf2, NULL, 0, 3);
 +}
 +
 +      i++;
 +}
 +
 +      CDDM_calc_edges(result);
 +      CDDM_calc_normals(result);
 +
 +      return result;
 +      */
 +      
 +      return NULL;
 +}
 +
 +
 +DerivedMesh *CDDM_create_tearing(ClothModifierData *clmd, DerivedMesh *dm)
 +{
 +      /*
 +      DerivedMesh *result = NULL;
 +      unsigned int i = 0, a = 0, j=0;
 +      int numverts = dm->getNumVerts(dm);
 +      int numedges = dm->getNumEdges(dm);
 +      int numfaces = dm->getNumFaces(dm);
 +
 +      MVert *mvert = CDDM_get_verts(dm);
 +      MEdge *medge = CDDM_get_edges(dm);
 +      MFace *mface = CDDM_get_faces(dm);
 +
 +      MVert *mvert2;
 +      MFace *mface2;
 +      unsigned int numtris=0;
 +      unsigned int numquads=0;
 +      EdgeHash *edgehash = NULL;
 +      Cloth *cloth = clmd->clothObject;
 +      ClothSpring *springs = cloth->springs;
 +      unsigned int numsprings = cloth->numsprings;
 +      
 +      // create spring tearing hash
 +      edgehash = BLI_edgehash_new();
 +      
 +      for(i = 0; i < numsprings; i++)
 +      {
 +      if((springs[i].flags & CSPRING_FLAG_DEACTIVATE)
 +      &&(!BLI_edgehash_haskey(edgehash, springs[i].ij, springs[i].kl)))
 +      {
 +      BLI_edgehash_insert(edgehash, springs[i].ij, springs[i].kl, NULL);
 +      BLI_edgehash_insert(edgehash, springs[i].kl, springs[i].ij, NULL);
 +      j++;
 +}
 +}
 +      
 +      // printf("found %d tears\n", j);
 +      
 +      result = CDDM_from_template(dm, numverts, 0, numfaces);
 +
 +      if(!result)
 +      return NULL;
 +
 +      // do verts
 +      mvert2 = CDDM_get_verts(result);
 +      for(a=0; a<numverts; a++) 
 +      {
 +      MVert *inMV;
 +      MVert *mv = &mvert2[a];
 +
 +      inMV = &mvert[a];
 +
 +      DM_copy_vert_data(dm, result, a, a, 1);
 +      *mv = *inMV;
 +}
 +
 +
 +      // do faces
 +      mface2 = CDDM_get_faces(result);
 +      for(a=0, i=0; a<numfaces; a++) 
 +      {
 +      MFace *mf = &mface2[i];
 +      MFace *inMF;
 +      inMF = &mface[a];
 +
 +              
 +              // DM_copy_face_data(dm, result, a, i, 1);
 +
 +              // *mf = *inMF;
 +              
 +              
 +      if((!BLI_edgehash_haskey(edgehash, mface[a].v1, mface[a].v2))
 +      &&(!BLI_edgehash_haskey(edgehash, mface[a].v2, mface[a].v3))
 +      &&(!BLI_edgehash_haskey(edgehash, mface[a].v3, mface[a].v4))
 +      &&(!BLI_edgehash_haskey(edgehash, mface[a].v4, mface[a].v1)))
 +      {
 +      mf->v1 = mface[a].v1;
 +      mf->v2 = mface[a].v2;
 +      mf->v3 = mface[a].v3;
 +      mf->v4 = mface[a].v4;
 +      
 +      test_index_face(mf, NULL, 0, 4);
 +      
 +      i++;
 +}
 +}
 +
 +      CDDM_lower_num_faces(result, i);
 +      CDDM_calc_edges(result);
 +      CDDM_calc_normals(result);
 +      
 +      BLI_edgehash_free(edgehash, NULL);
 +
 +      return result;
 +      */
 +      
 +      return NULL;
 +}
 +
 +int modifiers_indexInObject(Object *ob, ModifierData *md_seek);
 +
 +void cloth_clear_cache(Object *ob, ClothModifierData *clmd, float framenr)
 +{
 +      int stack_index = -1;
 +      
 +      if(!(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_CCACHE_PROTECT))
 +      {
 +              stack_index = modifiers_indexInObject(ob, (ModifierData *)clmd);
 +              
 +              BKE_ptcache_id_clear((ID *)ob, PTCACHE_CLEAR_AFTER, framenr, stack_index);
 +      }
 +}
 +static void cloth_write_cache(Object *ob, ClothModifierData *clmd, float framenr)
 +{
 +      FILE *fp = NULL;
 +      int stack_index = -1;
 +      unsigned int a;
 +      Cloth *cloth = clmd->clothObject;
 +      
 +      if(!cloth)
 +              return;
 +      
 +      stack_index = modifiers_indexInObject(ob, (ModifierData *)clmd);
 +      
 +      fp = BKE_ptcache_id_fopen((ID *)ob, 'w', framenr, stack_index);
 +      if(!fp) return;
 +      
 +      for(a = 0; a < cloth->numverts; a++)
 +      {
 +              fwrite(&cloth->x[a], sizeof(float),4,fp);
 +              fwrite(&cloth->xconst[a], sizeof(float),4,fp);
 +              fwrite(&cloth->v[a], sizeof(float),4,fp);
 +      }
 +      
 +      fclose(fp);
 +}
 +static int cloth_read_cache(Object *ob, ClothModifierData *clmd, float framenr)
 +{
 +      FILE *fp = NULL;
 +      int stack_index = -1;
 +      unsigned int a, ret = 1;
 +      Cloth *cloth = clmd->clothObject;
 +      
 +      if(!cloth)
 +              return 0;
 +      
 +      stack_index = modifiers_indexInObject(ob, (ModifierData *)clmd);
 +      
 +      fp = BKE_ptcache_id_fopen((ID *)ob, 'r', framenr, stack_index);
 +      if(!fp)
 +              ret = 0;
 +      else {
 +              for(a = 0; a < cloth->numverts; a++)
 +              {
 +                      if(fread(&cloth->x[a], sizeof(float), 4, fp) != 4) 
 +                      {
 +                              ret = 0;
 +                              break;
 +                      }
 +                      if(fread(&cloth->xconst[a], sizeof(float), 4, fp) != 4) 
 +                      {
 +                              ret = 0;
 +                              break;
 +                      }
 +                      if(fread(&cloth->v[a], sizeof(float), 4, fp) != 4) 
 +                      {
 +                              ret = 0;
 +                              break;
 +                      }
 +              }
 +              
 +              fclose(fp);
 +      }
 +      
 +      implicit_set_positions(clmd);
 +                      
 +      return ret;
 +}
 +
 +/************************************************
 + * clothModifier_do - main simulation function
 +************************************************/
 +DerivedMesh *clothModifier_do(ClothModifierData *clmd,Object *ob, DerivedMesh *dm, int useRenderParams, int isFinalCalc)
 +{
 +      unsigned int i;
 +      unsigned int numverts = -1;
 +      unsigned int numedges = -1;
 +      unsigned int numfaces = -1;
 +      MVert *mvert = NULL;
 +      MEdge *medge = NULL;
 +      MFace *mface = NULL;
 +      DerivedMesh *result = NULL;
 +      Cloth *cloth = clmd->clothObject;
 +      unsigned int framenr = (float)G.scene->r.cfra;
 +      float current_time = bsystem_time(ob, (float)G.scene->r.cfra, 0.0);
 +      ListBase *effectors = NULL;
 +      float deltaTime = current_time - clmd->sim_parms->sim_time;     
 +
 +      clmd->sim_parms->dt = 1.0f / (clmd->sim_parms->stepsPerFrame * G.scene->r.frs_sec);
 +
 +      result = CDDM_copy(dm);
 +      
 +      if(!result)
 +      {
 +              return dm;
 +      }
 +      
 +      numverts = result->getNumVerts(result);
 +      numedges = result->getNumEdges(result);
 +      numfaces = result->getNumFaces(result);
 +      mvert = CDDM_get_verts(result);
 +      medge = CDDM_get_edges(result);
 +      mface = CDDM_get_faces(result);
 +
 +      clmd->sim_parms->sim_time = current_time;
 +      
 +      if ( current_time < clmd->sim_parms->firstframe )
 +              return result;
 +      
 +      // only be active during a specific period:
 +      // that's "first frame" and "last frame" on GUI
 +      /*
 +      if ( clmd->clothObject )
 +      {
 +      if ( clmd->sim_parms->cache )
 +      {
 +      if ( current_time < clmd->sim_parms->firstframe )
 +      {
 +      int frametime = cloth_cache_first_frame ( clmd );
 +      if ( cloth_cache_search_frame ( clmd, frametime ) )
 +      {
 +      cloth_cache_get_frame ( clmd, frametime );
 +      cloth_to_object ( ob, result, clmd );
 +}
 +      return result;
 +}
 +      else if ( current_time > clmd->sim_parms->lastframe )
 +      {
 +      int frametime = cloth_cache_last_frame ( clmd );
 +      if ( cloth_cache_search_frame ( clmd, frametime ) )
 +      {
 +      cloth_cache_get_frame ( clmd, frametime );
 +      cloth_to_object ( ob, result, clmd );
 +}
 +      return result;
 +}
 +      else if ( ABS ( deltaTime ) >= 2.0f ) // no timewarps allowed
 +      {
 +      if ( cloth_cache_search_frame ( clmd, framenr ) )
 +      {
 +      cloth_cache_get_frame ( clmd, framenr );
 +      cloth_to_object ( ob, result, clmd );
 +}
 +      clmd->sim_parms->sim_time = current_time;
 +      return result;
 +}
 +}
 +}
 +      */
 +      
 +      if(deltaTime == 1.0f)
 +      {
 +              if ((clmd->clothObject == NULL) || (numverts != clmd->clothObject->numverts) ) 
 +              {
 +                      cloth_clear_cache(ob, clmd, 0);
 +                      
 +                      if(!cloth_from_object (ob, clmd, result, dm, framenr))
 +                              return result;
 +
 +                      if(clmd->clothObject == NULL)
 +                              return result;
 +
 +                      cloth = clmd->clothObject;
 +              }
 +
 +              clmd->clothObject->old_solver_type = clmd->sim_parms->solver_type;
 +
 +              // Insure we have a clmd->clothObject, in case allocation failed.
 +              if (clmd->clothObject != NULL) 
 +              {
 +                      if(!cloth_read_cache(ob, clmd, framenr))
 +                      {
 +                              // Force any pinned verts to their constrained location.
 +                              // has to be commented for verlet
 +                              for ( i = 0; i < clmd->clothObject->numverts; i++ )
 +                              {
 +                                      // Save the previous position.
 +                                      VECCOPY ( cloth->xold[i], cloth->xconst[i] );
 +                                      VECCOPY ( cloth->current_xold[i], cloth->x[i] );
 +                                      // Get the current position.
 +                                      VECCOPY ( cloth->xconst[i], mvert[i].co );
 +                                      Mat4MulVecfl ( ob->obmat, cloth->xconst[i] );
 +                              }
 +                              
 +                              tstart();
 +                              
 +                              /* Call the solver. */
 +                              if (solvers [clmd->sim_parms->solver_type].solver)
 +                                      solvers [clmd->sim_parms->solver_type].solver (ob, framenr, clmd, effectors);
 +                              
 +                              tend();
 +                              
 +                              printf("Cloth simulation time: %f\n", tval());
 +                              
 +                              cloth_write_cache(ob, clmd, framenr);
 +                      }
 +
 +                      // Copy the result back to the object.
 +                      cloth_to_object (ob, result, clmd);
 +                      
 +                      // bvh_free(clmd->clothObject->tree);
 +                      // clmd->clothObject->tree = bvh_build(clmd, clmd->coll_parms->epsilon);
 +              } 
 +
 +      }
 +      else if ( ( deltaTime <= 0.0f ) || ( deltaTime > 1.0f ) )
 +      {
 +              if ( clmd->clothObject != NULL )
 +              {
 +                      if(cloth_read_cache(ob, clmd, framenr))
 +                              cloth_to_object (ob, result, clmd);
 +              }
 +      }
 +      
 +      return result;
 +}
 +
 +/* frees all */
 +void cloth_free_modifier (ClothModifierData *clmd)
 +{
 +      Cloth   *cloth = NULL;
 +
 +      if(!clmd)
 +              return;
 +
 +      cloth = clmd->clothObject;
 +
 +      // free our frame cache
 +      // cloth_clear_cache(ob, clmd, 0);
 +
 +      /* Calls the solver and collision frees first as they
 +      * might depend on data in clmd->clothObject. */
 +
 +      if (cloth) 
 +      {       
 +              // If our solver provides a free function, call it
 +              if (cloth->old_solver_type < 255 && solvers [cloth->old_solver_type].free) 
 +              {       
 +                      solvers [cloth->old_solver_type].free (clmd);
 +              }
 +              
 +              // Free the verts.
 +              if (cloth->verts != NULL)
 +                      MEM_freeN (cloth->verts);
 +              
 +              // Free the faces.
 +              if ( cloth->mfaces != NULL )
 +                      MEM_freeN ( cloth->mfaces );
 +              
 +              // Free the verts.
 +              if ( cloth->x != NULL )
 +                      MEM_freeN ( cloth->x );
 +              
 +              // Free the verts.
 +              if ( cloth->xold != NULL )
 +                      MEM_freeN ( cloth->xold );
 +              
 +              // Free the verts.
 +              if ( cloth->v != NULL )
 +                      MEM_freeN ( cloth->v );
 +              
 +              // Free the verts.
 +              if ( cloth->current_x != NULL )
 +                      MEM_freeN ( cloth->current_x );
 +              
 +              // Free the verts.
 +              if ( cloth->current_xold != NULL )
 +                      MEM_freeN ( cloth->current_xold );
 +              
 +              // Free the verts.
 +              if ( cloth->current_v != NULL )
 +                      MEM_freeN ( cloth->current_v );
 +              
 +              // Free the verts.
 +              if ( cloth->xconst != NULL )
 +                      MEM_freeN ( cloth->xconst );
 +              
 +              cloth->verts = NULL;
 +              cloth->numverts = -1;
 +              
 +              // Free the springs.
 +              if ( cloth->springs != NULL )
 +              {
 +                      LinkNode *search = cloth->springs;
 +                      while(search)
 +                      {
 +                              ClothSpring *spring = search->link;
 +                                              
 +                              MEM_freeN ( spring );
 +                              search = search->next;
 +                      }
 +                      BLI_linklist_free(cloth->springs, NULL);
 +              
 +                      cloth->springs = NULL;
 +              }
 +
 +              cloth->numsprings = -1;
 +              
 +              // free BVH collision tree
 +              if(cloth->tree)
 +                      bvh_free((BVH *)cloth->tree);
 +              
 +              // free BVH self collision tree
 +              if(cloth->selftree)
 +                      bvh_free((BVH *)cloth->selftree);
 +              
 +              if(cloth->edgehash)
 +                      BLI_edgehash_free ( cloth->edgehash, NULL );
 +              
 +              MEM_freeN (cloth);
 +              clmd->clothObject = NULL;
 +      }
 +
 +}
 +
 +
 +
 +/******************************************************************************
 +*
 +* Internal functions.
 +*
 +******************************************************************************/
 +
 +/**
 + * cloth_to_object - copies the deformed vertices to the object.
 + *
 + * This function is a modified version of the softbody.c:softbody_to_object() function.
 + **/
 +static void cloth_to_object (Object *ob,  DerivedMesh *dm, ClothModifierData *clmd)
 +{
 +      unsigned int    i = 0;
 +      MVert *mvert = NULL;
 +      unsigned int numverts;
 +      Cloth *cloth = clmd->clothObject;
 +
 +      if (clmd->clothObject) {
 +              /* inverse matrix is not uptodate... */
 +              Mat4Invert (ob->imat, ob->obmat);
 +
 +              mvert = CDDM_get_verts(dm);
 +              numverts = dm->getNumVerts(dm);
 +
 +              for (i = 0; i < numverts; i++)
 +              {
 +                      VECCOPY (mvert[i].co, cloth->x[i]);
 +                      Mat4MulVecfl (ob->imat, mvert[i].co);   /* cloth is in global coords */
 +              }
 +      }
 +}
 +
 +
 +/**
 + * cloth_apply_vgroup - applies a vertex group as specified by type
 + *
 + **/
 +static void cloth_apply_vgroup(ClothModifierData *clmd, DerivedMesh *dm, short vgroup)
 +{
 +      unsigned int i = 0;
 +      unsigned int j = 0;
 +      MDeformVert     *dvert = NULL;
 +      Cloth *clothObj = NULL;
 +      unsigned int numverts = 0;
 +      float goalfac = 0;
 +      ClothVertex *verts = NULL;
 +
 +      clothObj = clmd->clothObject;
 +
 +      if ( !dm )
 +              return;
 +
 +      numverts = dm->getNumVerts(dm);
 +
 +      /* vgroup is 1 based, decrement so we can match the right group. */
 +      --vgroup;
 +
 +      verts = clothObj->verts;
 +
 +      for ( i = 0; i < numverts; i++, verts++ )
 +      {
 +              // LATER ON, support also mass painting here
 +              if ( clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL )
 +              {
 +                      dvert = dm->getVertData ( dm, i, CD_MDEFORMVERT );
 +                      if ( dvert )
 +                      {
 +                              for ( j = 0; j < dvert->totweight; j++ )
 +                              {
 +                                      if ( dvert->dw[j].def_nr == vgroup )
 +                                      {
 +                                              verts->goal = dvert->dw [j].weight;
 +
 +                                              goalfac= ABS ( clmd->sim_parms->maxgoal - clmd->sim_parms->mingoal );
 +                                              verts->goal  = ( float ) pow ( verts->goal , 4.0f );
 +
 +                                              if ( dvert->dw [j].weight >=SOFTGOALSNAP )
 +                                              {
 +                                                      verts->flags |= CVERT_FLAG_PINNED;
 +                                              }
 +
 +                                              // TODO enable mass painting here, for the moment i let "goals" go first
 +
 +                                              break;
 +                                      }
 +                              }
 +                      }
 +              }
 +      }
 +}
 +              
 +// only meshes supported at the moment
 +static int cloth_from_object(Object *ob, ClothModifierData *clmd, DerivedMesh *dm, DerivedMesh *olddm, float framenr)
 +{
 +      unsigned int i;
 +      unsigned int numverts = dm->getNumVerts(dm);
 +      MVert *mvert = CDDM_get_verts(dm);
 +      float tnull[3] = {0,0,0};
 +      
 +      /* If we have a clothObject, free it. */
 +      if (clmd->clothObject != NULL)
 +              cloth_free_modifier (clmd);
 +
 +      /* Allocate a new cloth object. */
 +      clmd->clothObject = MEM_callocN (sizeof(Cloth), "cloth");
 +      if (clmd->clothObject) 
 +      {
 +              clmd->clothObject->old_solver_type = 255;
 +              clmd->clothObject->edgehash = NULL;
 +      }
 +      else if (clmd->clothObject == NULL) 
 +      {
 +              modifier_setError (&(clmd->modifier), "Out of memory on allocating clmd->clothObject.");
 +              return 0;
 +      }
 +
 +      switch (ob->type)
 +      {
 +              case OB_MESH:
 +              
 +                      // mesh input objects need DerivedMesh
 +                      if ( !dm )
 +                              return 0;
 +
 +                      cloth_from_mesh (ob, clmd, dm, framenr);
 +
 +                      if ( clmd->clothObject != NULL )
 +                      {
 +                              /* create springs */
 +                              clmd->clothObject->springs = NULL;
 +                              clmd->clothObject->numsprings = -1;
 +                                      
 +                              /* set initial values */
 +                              for (i = 0; i < numverts; ++i)
 +                              {
 +                                      VECCOPY (clmd->clothObject->x[i], mvert[i].co);
 +                                      Mat4MulVecfl(ob->obmat, clmd->clothObject->x[i]);
 +      
 +                                      clmd->clothObject->verts [i].mass = clmd->sim_parms->mass;
 +                                      if ( clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL )
 +                                              clmd->clothObject->verts [i].goal= clmd->sim_parms->defgoal;
 +                                      else
 +                                              clmd->clothObject->verts [i].goal= 0.0;
 +                                      clmd->clothObject->verts [i].flags = 0;
 +                                      VECCOPY(clmd->clothObject->xold[i], clmd->clothObject->x[i]);
 +                                      VECCOPY(clmd->clothObject->xconst[i], clmd->clothObject->x[i]);
 +                                      VECCOPY(clmd->clothObject->current_xold[i], clmd->clothObject->x[i]);
 +                                      VecMulf(clmd->clothObject->v[i], 0.0);
 +      
 +                                      clmd->clothObject->verts [i].impulse_count = 0;
 +                                      VECCOPY ( clmd->clothObject->verts [i].impulse, tnull );
 +                              }
 +                              
- int cloth_build_springs ( Cloth *cloth, DerivedMesh *dm )
++                              if (!cloth_build_springs (clmd, dm) )
 +                              {
 +                                      modifier_setError (&(clmd->modifier), "Can't build springs.");
 +                                      return 0;
 +                              }  
 +      
 +                              /* apply / set vertex groups */
 +                              if (clmd->sim_parms->vgroup_mass > 0)
 +                                      cloth_apply_vgroup (clmd, olddm, clmd->sim_parms->vgroup_mass);
 +      
 +                              /* init our solver */
 +                              if (solvers [clmd->sim_parms->solver_type].init)
 +                                      solvers [clmd->sim_parms->solver_type].init (ob, clmd);
 +      
 +                              clmd->clothObject->tree = bvh_build_from_float3(CDDM_get_faces(dm), dm->getNumFaces(dm), clmd->clothObject->x, numverts, clmd->coll_parms->epsilon);
 +                              
 +                              clmd->clothObject->selftree = bvh_build_from_float3(NULL, 0, clmd->clothObject->x, numverts, clmd->coll_parms->selfepsilon);
 +                              
 +                              // save initial state
 +                              cloth_write_cache(ob, clmd, framenr-1);
 +                      }
 +                      return 1;
 +                      default: return 0; // TODO - we do not support changing meshes
 +      }
 +      
 +      return 0;
 +}
 +
 +static void cloth_from_mesh (Object *ob, ClothModifierData *clmd, DerivedMesh *dm, float framenr)
 +{
 +      unsigned int numverts = dm->getNumVerts(dm);
 +      unsigned int numfaces = dm->getNumFaces(dm);
 +      MFace *mface = CDDM_get_faces(dm);
 +
 +      /* Allocate our vertices.
 +      */
 +      clmd->clothObject->numverts = numverts;
 +      clmd->clothObject->verts = MEM_callocN (sizeof (ClothVertex) * clmd->clothObject->numverts, "clothVertex");
 +      if (clmd->clothObject->verts == NULL) 
 +      {
 +              cloth_free_modifier (clmd);
 +              modifier_setError (&(clmd->modifier), "Out of memory on allocating clmd->clothObject->verts.");
 +              return;
 +      }
 +      
 +      clmd->clothObject->x = MEM_callocN ( sizeof ( float ) * clmd->clothObject->numverts * 4, "Cloth MVert_x" );
 +      if ( clmd->clothObject->x == NULL )
 +      {
 +              cloth_free_modifier ( clmd );
 +              modifier_setError ( & ( clmd->modifier ), "Out of memory on allocating clmd->clothObject->x." );
 +              return;
 +      }
 +      
 +      clmd->clothObject->xold = MEM_callocN ( sizeof ( float ) * clmd->clothObject->numverts * 4, "Cloth MVert_xold" );
 +      if ( clmd->clothObject->xold == NULL )
 +      {
 +              cloth_free_modifier ( clmd );
 +              modifier_setError ( & ( clmd->modifier ), "Out of memory on allocating clmd->clothObject->xold." );
 +              return;
 +      }
 +      
 +      clmd->clothObject->v = MEM_callocN ( sizeof ( float ) * clmd->clothObject->numverts * 4, "Cloth MVert_v" );
 +      if ( clmd->clothObject->v == NULL )
 +      {
 +              cloth_free_modifier ( clmd );
 +              modifier_setError ( & ( clmd->modifier ), "Out of memory on allocating clmd->clothObject->v." );
 +              return;
 +      }
 +      
 +      clmd->clothObject->current_x = MEM_callocN ( sizeof ( float ) * clmd->clothObject->numverts * 4, "Cloth MVert_current_x" );
 +      if ( clmd->clothObject->current_x == NULL )
 +      {
 +              cloth_free_modifier ( clmd );
 +              modifier_setError ( & ( clmd->modifier ), "Out of memory on allocating clmd->clothObject->current_x." );
 +              return;
 +      }
 +      
 +      clmd->clothObject->current_xold = MEM_callocN ( sizeof ( float ) * clmd->clothObject->numverts * 4, "Cloth MVert_current_xold" );
 +      if ( clmd->clothObject->current_xold == NULL )
 +      {
 +              cloth_free_modifier ( clmd );
 +              modifier_setError ( & ( clmd->modifier ), "Out of memory on allocating clmd->clothObject->current_xold." );
 +              return;
 +      }
 +      
 +      clmd->clothObject->current_v = MEM_callocN ( sizeof ( float ) * clmd->clothObject->numverts * 4, "Cloth MVert_current_v" );
 +      if ( clmd->clothObject->current_v == NULL )
 +      {
 +              cloth_free_modifier ( clmd );
 +              modifier_setError ( & ( clmd->modifier ), "Out of memory on allocating clmd->clothObject->current_v." );
 +              return;
 +      }
 +      
 +      clmd->clothObject->xconst = MEM_callocN ( sizeof ( float ) * clmd->clothObject->numverts * 4, "Cloth MVert_xconst" );
 +      if ( clmd->clothObject->xconst == NULL )
 +      {
 +              cloth_free_modifier ( clmd );
 +              modifier_setError ( & ( clmd->modifier ), "Out of memory on allocating clmd->clothObject->xconst." );
 +              return;
 +      }
 +
 +      // save face information
 +      clmd->clothObject->numfaces = numfaces;
 +      clmd->clothObject->mfaces = MEM_callocN (sizeof (MFace) * clmd->clothObject->numfaces, "clothMFaces");
 +      if (clmd->clothObject->mfaces == NULL) 
 +      {
 +              cloth_free_modifier (clmd);
 +              modifier_setError (&(clmd->modifier), "Out of memory on allocating clmd->clothObject->mfaces.");
 +              return;
 +      }
 +      memcpy(clmd->clothObject->mfaces, mface, sizeof(MFace)*numfaces);
 +
 +      /* Free the springs since they can't be correct if the vertices
 +      * changed.
 +      */
 +      if (clmd->clothObject->springs != NULL)
 +              MEM_freeN (clmd->clothObject->springs);
 +
 +}
 +
 +/***************************************************************************************
 +* SPRING NETWORK BUILDING IMPLEMENTATION BEGIN
 +***************************************************************************************/
 +
 +// be carefull: implicit solver has to be resettet when using this one!
 +int cloth_add_spring ( ClothModifierData *clmd, unsigned int indexA, unsigned int indexB, float restlength, int spring_type)
 +{
 +      Cloth *cloth = clmd->clothObject;
 +      ClothSpring *spring = NULL;
 +      
 +      if(cloth)
 +      {
 +              // TODO: look if this spring is already there
 +              
 +              spring = ( ClothSpring * ) MEM_callocN ( sizeof ( ClothSpring ), "cloth spring" );
 +              
 +              spring->ij = indexA;
 +              spring->kl = indexB;
 +              spring->restlen =  restlength;
 +              spring->type = spring_type;
 +              spring->flags = 0;
 +              
 +              cloth->numsprings++;
 +      
 +              BLI_linklist_append ( &cloth->springs, spring );
 +              
 +              return 1;
 +      }
 +      return 0;
 +}
 +
-       unsigned int i = 0, j = 0;
++int cloth_build_springs ( ClothModifierData *clmd, DerivedMesh *dm )
 +{
++      Cloth *cloth = clmd->clothObject;
 +      ClothSpring *spring = NULL, *tspring = NULL, *tspring2 = NULL;
 +      unsigned int struct_springs = 0, shear_springs=0, bend_springs = 0;
-       float temp[3];
++      unsigned int i = 0, j = 0, akku_count;
 +      unsigned int numverts = dm->getNumVerts ( dm );
 +      unsigned int numedges = dm->getNumEdges ( dm );
 +      unsigned int numfaces = dm->getNumFaces ( dm );
 +      MEdge *medge = CDDM_get_edges ( dm );
 +      MFace *mface = CDDM_get_faces ( dm );
 +      unsigned int index2 = 0; // our second vertex index
 +      LinkNode **edgelist = NULL;
 +      EdgeHash *edgehash = NULL;
 +      LinkNode *search = NULL, *search2 = NULL;
-       /*
++      float temp[3], akku, min, max;
 +      LinkNode *node = NULL, *node2 = NULL;
 +      
 +      // error handling
 +      if ( numedges==0 )
 +              return 0;
 +
 +      cloth->springs = NULL;
 +
 +      edgelist = MEM_callocN ( sizeof ( LinkNode * ) * numverts, "cloth_edgelist_alloc" );
 +      for ( i = 0; i < numverts; i++ )
 +      {
 +              edgelist[i] = NULL;
 +      }
 +
 +      if ( cloth->springs )
 +              MEM_freeN ( cloth->springs );
 +
 +      // create spring network hash
 +      edgehash = BLI_edgehash_new();
 +
 +      // structural springs
 +      for ( i = 0; i < numedges; i++ )
 +      {
 +              spring = ( ClothSpring * ) MEM_callocN ( sizeof ( ClothSpring ), "cloth spring" );
 +
 +              if ( spring )
 +              {
 +                      spring->ij = medge[i].v1;
 +                      spring->kl = medge[i].v2;
 +                      VECSUB ( temp, cloth->x[spring->kl], cloth->x[spring->ij] );
 +                      spring->restlen =  sqrt ( INPR ( temp, temp ) );
 +                      spring->type = CLOTH_SPRING_TYPE_STRUCTURAL;
 +                      spring->flags = 0;
 +                      struct_springs++;
 +                      
 +                      if(!i)
 +                              node2 = BLI_linklist_append_fast ( &cloth->springs, spring );
 +                      else
 +                              node2 = BLI_linklist_append_fast ( &node->next, spring );
 +                      node = node2;
 +              }
 +      }
 +      
 +      // calc collision balls *slow*
 +      // better: use precalculated list with O(1) index access to all springs of a vertex
 +      // missing for structural since it's not needed for building bending springs
-                               akku += bs->len;
-                               akku_count++,
-                               min = MIN2(bs->len,min);
-                               max = MAX2(bs->len,max);
 +      for ( i = 0; i < numverts; i++ )
 +      {
++              akku_count = 0;
++              akku = 0.0;
++              cloth->verts[i].collball=0;
++              min = 1e22f;
++              max = -1e22f;
++              
 +              search = cloth->springs;
 +              for ( j = 0; j < struct_springs; j++ )
 +              {
 +                      if ( !search )
 +                              break;
 +
 +                      tspring = search->link;
 +                      
 +                      if((tspring->ij == i) || (tspring->kl == i))
 +                      {
-       */
++                              akku += spring->restlen;
++                              akku_count++;
++                              min = MIN2(spring->restlen,min);
++                              max = MAX2(spring->restlen,max);
 +                      }
 +              }
++              
++              if (akku_count > 0) {
++                      cloth->verts[i].collball = akku/(float)akku_count*clmd->coll_parms->selfepsilon;
++              }
++              else cloth->verts[i].collball=0;
 +      }
++      
 +      // shear springs
 +      for ( i = 0; i < numfaces; i++ )
 +      {
 +              spring = ( ClothSpring *) MEM_callocN ( sizeof ( ClothSpring ), "cloth spring" );
 +
 +              spring->ij = mface[i].v1;
 +              spring->kl = mface[i].v3;
 +              VECSUB ( temp, cloth->x[spring->kl], cloth->x[spring->ij] );
 +              spring->restlen =  sqrt ( INPR ( temp, temp ) );
 +              spring->type = CLOTH_SPRING_TYPE_SHEAR;
 +
 +              BLI_linklist_append ( &edgelist[spring->ij], spring );
 +              BLI_linklist_append ( &edgelist[spring->kl], spring );
 +              shear_springs++;
 +
 +              node2 = BLI_linklist_append_fast ( &node->next, spring );
 +              node = node2;
 +
 +              if ( mface[i].v4 )
 +              {
 +                      spring = ( ClothSpring * ) MEM_callocN ( sizeof ( ClothSpring ), "cloth spring" );
 +
 +                      spring->ij = mface[i].v2;
 +                      spring->kl = mface[i].v4;
 +                      VECSUB ( temp, cloth->x[spring->kl], cloth->x[spring->ij] );
 +                      spring->restlen =  sqrt ( INPR ( temp, temp ) );
 +                      spring->type = CLOTH_SPRING_TYPE_SHEAR;
 +
 +                      BLI_linklist_append ( &edgelist[spring->ij], spring );
 +                      BLI_linklist_append ( &edgelist[spring->kl], spring );
 +                      shear_springs++;
 +
 +                      node2 = BLI_linklist_append_fast ( &node->next, spring );
 +                      node = node2;
 +              }
 +      }
 +      
 +      // bending springs
 +      search2 = cloth->springs;
 +      for ( i = struct_springs; i < struct_springs+shear_springs; i++ )
 +      {
 +              if ( !search2 )
 +                      break;
 +
 +              tspring2 = search2->link;
 +              search = edgelist[tspring2->kl];
 +              while ( search )
 +              {
 +                      tspring = search->link;
 +                      index2 = ( ( tspring->ij==tspring2->kl ) ? ( tspring->kl ) : ( tspring->ij ) );
 +                      
 +                      // check for existing spring
 +                      // check also if startpoint is equal to endpoint
 +                      if ( !BLI_edgehash_haskey ( edgehash, index2, tspring2->ij )
 +                      && !BLI_edgehash_haskey ( edgehash, tspring2->ij, index2 )
 +                      && ( index2!=tspring2->ij ) )
 +                      {
 +                              spring = ( ClothSpring * ) MEM_callocN ( sizeof ( ClothSpring ), "cloth spring" );
 +
 +                              spring->ij = tspring2->ij;
 +                              spring->kl = index2;
 +                              VECSUB ( temp, cloth->x[index2], cloth->x[tspring2->ij] );
 +                              spring->restlen =  sqrt ( INPR ( temp, temp ) );
 +                              spring->type = CLOTH_SPRING_TYPE_BENDING;
 +                              BLI_edgehash_insert ( edgehash, spring->ij, index2, NULL );
 +                              bend_springs++;
 +
 +                              node2 = BLI_linklist_append_fast ( &node->next, spring );
 +                              node = node2;
 +                      }
 +                      search = search->next;
 +              }
 +              search2 = search2->next;
 +      }
 +      
 +      cloth->numspringssave = cloth->numsprings = struct_springs + shear_springs + bend_springs;
 +      cloth->numothersprings = struct_springs + shear_springs;
 +      
 +      for ( i = 0; i < numverts; i++ )
 +      {
 +              BLI_linklist_free ( edgelist[i],NULL );
 +      }
 +      if ( edgelist )
 +              MEM_freeN ( edgelist );
 +      
 +      cloth->edgehash = edgehash;
 +
 +      return 1;
 +
 +} /* cloth_build_springs */
 +/***************************************************************************************
 +* SPRING NETWORK BUILDING IMPLEMENTATION END
 +***************************************************************************************/
 +
index 0f8357cee251c74a3ef5318eee6dc53204b6d876,0000000000000000000000000000000000000000..12c5c1906911c262ae9d5f8029a65280fd413057
mode 100644,000000..100644
--- /dev/null
@@@ -1,2310 -1,0 +1,2315 @@@
- /*
- #pragma omp parallel for private(i,j, collisions) shared(current_x)
-       for(count = 0; count < 6; count++)
-       {
-       collisions = 0;
-               
-       for(i = 0; i < cloth->numverts; i++)
-       {
-       for(j = i + 1; j < cloth->numverts; j++)
-       {
-       float temp[3];
-       float length = 0;
-       float mindistance = cloth->selftree->epsilon;
-                                       
-       if(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL)
-       {                       
-       if((cloth->verts [i].goal >= SOFTGOALSNAP)
-       && (cloth->verts [j].goal >= SOFTGOALSNAP))
-       {
-       continue;
- }
- }
-                                       
-                                       // check for adjacent points
-       if(BLI_edgehash_haskey ( cloth->edgehash, i, j ))
-       {
-       continue;
- }
-                                       
-       VECSUB(temp, current_x[i], current_x[j]);
-                                       
-       length = Normalize(temp);
-                                       
-       if(length < mindistance)
-       {
-       float correction = mindistance - length;
-                                               
-       if((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [i].goal >= SOFTGOALSNAP))
-       {
-       VecMulf(temp, -correction);
-       VECADD(current_x[j], current_x[j], temp);
- }
-       else if((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [j].goal >= SOFTGOALSNAP))
-       {
-       VecMulf(temp, correction);
-       VECADD(current_x[i], current_x[i], temp);
- }
-       else
 +/*  implicit.c      
 +* 
 +*
 +* ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
 +*
 +* This program is free software; you can redistribute it and/or
 +* modify it under the terms of the GNU General Public License
 +* as published by the Free Software Foundation; either version 2
 +* of the License, or (at your option) any later version. The Blender
 +* Foundation also sells licenses for use in proprietary software under
 +* the Blender License.  See http://www.blender.org/BL/ for information
 +* about this.
 +*
 +* This program is distributed in the hope that it will be useful,
 +* but WITHOUT ANY WARRANTY; without even the implied warranty of
 +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 +* GNU General Public License for more details.
 +*
 +* You should have received a copy of the GNU General Public License
 +* along with this program; if not, write to the Free Software Foundation,
 +* Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
 +*
 +* The Original Code is Copyright (C) Blender Foundation
 +* All rights reserved.
 +*
 +* The Original Code is: all of this file.
 +*
 +* Contributor(s): none yet.
 +*
 +* ***** END GPL/BL DUAL LICENSE BLOCK *****
 +*/
 +#include <math.h>
 +#include <stdlib.h>
 +#include <string.h>
 +#include <stdio.h>
 +#include "MEM_guardedalloc.h"
 +/* types */
 +#include "DNA_curve_types.h"
 +#include "DNA_object_types.h"
 +#include "DNA_object_force.h" 
 +#include "DNA_cloth_types.h"  
 +#include "DNA_key_types.h"
 +#include "DNA_mesh_types.h"
 +#include "DNA_modifier_types.h"
 +#include "DNA_meshdata_types.h"
 +#include "DNA_lattice_types.h"
 +#include "DNA_scene_types.h"
 +#include "DNA_modifier_types.h"
 +#include "BLI_arithb.h"
 +#include "BLI_blenlib.h"
 +#include "BLI_edgehash.h"
 +#include "BLI_threads.h"
 +#include "BKE_collisions.h"
 +#include "BKE_curve.h"
 +#include "BKE_displist.h"
 +#include "BKE_effect.h"
 +#include "BKE_global.h"
 +#include "BKE_key.h"
 +#include "BKE_object.h"
 +#include "BKE_cloth.h"
 +#include "BKE_modifier.h"
 +#include "BKE_utildefines.h"
 +#include "BKE_global.h"
 +#include  "BIF_editdeform.h"
 +
 +#ifdef _WIN32
 +#include <windows.h>
 +static LARGE_INTEGER _itstart, _itend;
 +static LARGE_INTEGER ifreq;
 +void itstart(void)
 +{
 +      static int first = 1;
 +      if(first) {
 +              QueryPerformanceFrequency(&ifreq);
 +              first = 0;
 +      }
 +      QueryPerformanceCounter(&_itstart);
 +}
 +void itend(void)
 +{
 +      QueryPerformanceCounter(&_itend);
 +}
 +double itval()
 +{
 +      return ((double)_itend.QuadPart -
 +                      (double)_itstart.QuadPart)/((double)ifreq.QuadPart);
 +}
 +#else
 +#include <sys/time.h>
 +
 +                       static struct timeval _itstart, _itend;
 +       static struct timezone itz;
 +       void itstart(void)
 +{
 +      gettimeofday(&_itstart, &itz);
 +}
 +void itend(void)
 +{
 +      gettimeofday(&_itend,&itz);
 +}
 +double itval()
 +{
 +      double t1, t2;
 +      t1 =  (double)_itstart.tv_sec + (double)_itstart.tv_usec/(1000*1000);
 +      t2 =  (double)_itend.tv_sec + (double)_itend.tv_usec/(1000*1000);
 +      return t2-t1;
 +}
 +#endif
 +
 +struct Cloth;
 +
 +//////////////////////////////////////////
 +/* fast vector / matrix library, enhancements are welcome :) -dg */
 +/////////////////////////////////////////
 +
 +/* DEFINITIONS */
 +typedef float lfVector[3];
 +typedef struct fmatrix3x3 {
 +      float m[3][3]; /* 4x4 matrix */
 +      unsigned int c,r; /* column and row number */
 +      int pinned; /* is this vertex allowed to move? */
 +      float n1,n2,n3; /* three normal vectors for collision constrains */
 +      unsigned int vcount; /* vertex count */
 +      unsigned int scount; /* spring count */ 
 +} fmatrix3x3;
 +
 +///////////////////////////
 +// float[3] vector
 +///////////////////////////
 +/* simple vector code */
 +/* STATUS: verified */
 +DO_INLINE void mul_fvector_S(float to[3], float from[3], float scalar)
 +{
 +      to[0] = from[0] * scalar;
 +      to[1] = from[1] * scalar;
 +      to[2] = from[2] * scalar;
 +}
 +/* simple cross product */
 +/* STATUS: verified */
 +DO_INLINE void cross_fvector(float to[3], float vectorA[3], float vectorB[3])
 +{
 +      to[0] = vectorA[1] * vectorB[2] - vectorA[2] * vectorB[1];
 +      to[1] = vectorA[2] * vectorB[0] - vectorA[0] * vectorB[2];
 +      to[2] = vectorA[0] * vectorB[1] - vectorA[1] * vectorB[0];
 +}
 +/* simple v^T * v product ("outer product") */
 +/* STATUS: HAS TO BE verified (*should* work) */
 +DO_INLINE void mul_fvectorT_fvector(float to[3][3], float vectorA[3], float vectorB[3])
 +{
 +      mul_fvector_S(to[0], vectorB, vectorA[0]);
 +      mul_fvector_S(to[1], vectorB, vectorA[1]);
 +      mul_fvector_S(to[2], vectorB, vectorA[2]);
 +}
 +/* simple v^T * v product with scalar ("outer product") */
 +/* STATUS: HAS TO BE verified (*should* work) */
 +DO_INLINE void mul_fvectorT_fvectorS(float to[3][3], float vectorA[3], float vectorB[3], float aS)
 +{
 +      mul_fvector_S(to[0], vectorB, vectorA[0]* aS);
 +      mul_fvector_S(to[1], vectorB, vectorA[1]* aS);
 +      mul_fvector_S(to[2], vectorB, vectorA[2]* aS);
 +}
 +
 +/* printf vector[3] on console: for debug output */
 +void print_fvector(float m3[3])
 +{
 +      printf("%f\n%f\n%f\n\n",m3[0],m3[1],m3[2]);
 +}
 +
 +///////////////////////////
 +// long float vector float (*)[3]
 +///////////////////////////
 +/* print long vector on console: for debug output */
 +DO_INLINE void print_lfvector(lfVector *fLongVector, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +      for(i = 0; i < verts; i++)
 +      {
 +              print_fvector(fLongVector[i]);
 +      }
 +}
 +/* create long vector */
 +DO_INLINE lfVector *create_lfvector(unsigned int verts)
 +{
 +      // TODO: check if memory allocation was successfull */
 +      return  (lfVector *)MEM_callocN (verts * sizeof(lfVector), "cloth_implicit_alloc_vector");
 +      // return (lfVector *)cloth_aligned_malloc(&MEMORY_BASE, verts * sizeof(lfVector));
 +}
 +/* delete long vector */
 +DO_INLINE void del_lfvector(lfVector *fLongVector)
 +{
 +      if (fLongVector != NULL)
 +      {
 +              MEM_freeN (fLongVector);
 +              // cloth_aligned_free(&MEMORY_BASE, fLongVector);
 +      }
 +}
 +/* copy long vector */
 +DO_INLINE void cp_lfvector(float (*to)[3], float (*from)[3], unsigned int verts)
 +{
 +      memcpy(to, from, verts * sizeof(lfVector));
 +}
 +/* init long vector with float[3] */
 +DO_INLINE void init_lfvector(lfVector *fLongVector, float vector[3], unsigned int verts)
 +{
 +      unsigned int i = 0;
 +      for(i = 0; i < verts; i++)
 +      {
 +              VECCOPY(fLongVector[i], vector);
 +      }
 +}
 +/* zero long vector with float[3] */
 +DO_INLINE void zero_lfvector(float (*to)[3], unsigned int verts)
 +{
 +      memset(to, 0.0f, verts * sizeof(lfVector));
 +}
 +/* multiply long vector with scalar*/
 +DO_INLINE void mul_lfvectorS(float (*to)[3], lfVector *fLongVector, float scalar, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +
 +      for(i = 0; i < verts; i++)
 +      {
 +              mul_fvector_S(to[i], fLongVector[i], scalar);
 +      }
 +}
 +/* multiply long vector with scalar*/
 +/* A -= B * float */
 +DO_INLINE void submul_lfvectorS(float (*to)[3], lfVector *fLongVector, float scalar, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +      for(i = 0; i < verts; i++)
 +      {
 +              VECSUBMUL(to[i], fLongVector[i], scalar);
 +      }
 +}
 +/* dot product for big vector */
 +DO_INLINE float dot_lfvector(lfVector *fLongVectorA, lfVector *fLongVectorB, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +      float temp = 0.0;
 +// schedule(guided, 2)
 +#pragma omp parallel for reduction(+: temp)
 +      for(i = 0; i < verts; i++)
 +      {
 +              temp += INPR(fLongVectorA[i], fLongVectorB[i]);
 +      }
 +      return temp;
 +}
 +/* A = B + C  --> for big vector */
 +DO_INLINE void add_lfvector_lfvector(float (*to)[3], lfVector *fLongVectorA, lfVector *fLongVectorB, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +
 +      for(i = 0; i < verts; i++)
 +      {
 +              VECADD(to[i], fLongVectorA[i], fLongVectorB[i]);
 +      }
 +
 +}
 +/* A = B + C * float --> for big vector */
 +DO_INLINE void add_lfvector_lfvectorS(float (*to)[3], lfVector *fLongVectorA, lfVector *fLongVectorB, float bS, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +
 +      for(i = 0; i < verts; i++)
 +      {
 +              VECADDS(to[i], fLongVectorA[i], fLongVectorB[i], bS);
 +
 +      }
 +}
 +/* A = B * float + C * float --> for big vector */
 +DO_INLINE void add_lfvectorS_lfvectorS(float (*to)[3], lfVector *fLongVectorA, float aS, lfVector *fLongVectorB, float bS, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +
 +      for(i = 0; i < verts; i++)
 +      {
 +              VECADDSS(to[i], fLongVectorA[i], aS, fLongVectorB[i], bS);
 +      }
 +}
 +/* A = B - C * float --> for big vector */
 +DO_INLINE void sub_lfvector_lfvectorS(float (*to)[3], lfVector *fLongVectorA, lfVector *fLongVectorB, float bS, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +      for(i = 0; i < verts; i++)
 +      {
 +              VECSUBS(to[i], fLongVectorA[i], fLongVectorB[i], bS);
 +      }
 +
 +}
 +/* A = B - C --> for big vector */
 +DO_INLINE void sub_lfvector_lfvector(float (*to)[3], lfVector *fLongVectorA, lfVector *fLongVectorB, unsigned int verts)
 +{
 +      unsigned int i = 0;
 +
 +      for(i = 0; i < verts; i++)
 +      {
 +              VECSUB(to[i], fLongVectorA[i], fLongVectorB[i]);
 +      }
 +
 +}
 +///////////////////////////
 +// 4x4 matrix
 +///////////////////////////
 +/* printf 4x4 matrix on console: for debug output */
 +void print_fmatrix(float m3[3][3])
 +{
 +      printf("%f\t%f\t%f\n",m3[0][0],m3[0][1],m3[0][2]);
 +      printf("%f\t%f\t%f\n",m3[1][0],m3[1][1],m3[1][2]);
 +      printf("%f\t%f\t%f\n\n",m3[2][0],m3[2][1],m3[2][2]);
 +}
 +/* copy 4x4 matrix */
 +DO_INLINE void cp_fmatrix(float to[3][3], float from[3][3])
 +{
 +      // memcpy(to, from, sizeof (float) * 9);
 +      VECCOPY(to[0], from[0]);
 +      VECCOPY(to[1], from[1]);
 +      VECCOPY(to[2], from[2]);
 +}
 +/* calculate determinant of 4x4 matrix */
 +DO_INLINE float det_fmatrix(float m[3][3])
 +{
 +      return  m[0][0]*m[1][1]*m[2][2] + m[1][0]*m[2][1]*m[0][2] + m[0][1]*m[1][2]*m[2][0] 
 +                      -m[0][0]*m[1][2]*m[2][1] - m[0][1]*m[1][0]*m[2][2] - m[2][0]*m[1][1]*m[0][2];
 +}
 +DO_INLINE void inverse_fmatrix(float to[3][3], float from[3][3])
 +{
 +      unsigned int i, j;
 +      float d;
 +
 +      if((d=det_fmatrix(from))==0)
 +      {
 +              printf("can't build inverse");
 +              exit(0);
 +      }
 +      for(i=0;i<3;i++) 
 +      {
 +              for(j=0;j<3;j++) 
 +              {
 +                      int i1=(i+1)%3;
 +                      int i2=(i+2)%3;
 +                      int j1=(j+1)%3;
 +                      int j2=(j+2)%3;
 +                      // reverse indexs i&j to take transpose
 +                      to[j][i] = (from[i1][j1]*from[i2][j2]-from[i1][j2]*from[i2][j1])/d;
 +                      /*
 +                      if(i==j)
 +                      to[i][j] = 1.0f / from[i][j];
 +                      else
 +                      to[i][j] = 0;
 +                      */
 +              }
 +      }
 +
 +}
 +
 +/* 4x4 matrix multiplied by a scalar */
 +/* STATUS: verified */
 +DO_INLINE void mul_fmatrix_S(float matrix[3][3], float scalar)
 +{
 +      mul_fvector_S(matrix[0], matrix[0],scalar);
 +      mul_fvector_S(matrix[1], matrix[1],scalar);
 +      mul_fvector_S(matrix[2], matrix[2],scalar);
 +}
 +
 +/* a vector multiplied by a 4x4 matrix */
 +/* STATUS: verified */
 +DO_INLINE void mul_fvector_fmatrix(float *to, float *from, float matrix[3][3])
 +{
 +      to[0] = matrix[0][0]*from[0] + matrix[1][0]*from[1] + matrix[2][0]*from[2];
 +      to[1] = matrix[0][1]*from[0] + matrix[1][1]*from[1] + matrix[2][1]*from[2];
 +      to[2] = matrix[0][2]*from[0] + matrix[1][2]*from[1] + matrix[2][2]*from[2];
 +}
 +
 +/* 4x4 matrix multiplied by a vector */
 +/* STATUS: verified */
 +DO_INLINE void mul_fmatrix_fvector(float *to, float matrix[3][3], float *from)
 +{
 +      to[0] = INPR(matrix[0],from);
 +      to[1] = INPR(matrix[1],from);
 +      to[2] = INPR(matrix[2],from);
 +}
 +/* 4x4 matrix multiplied by a 4x4 matrix */
 +/* STATUS: verified */
 +DO_INLINE void mul_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      mul_fvector_fmatrix(to[0], matrixA[0],matrixB);
 +      mul_fvector_fmatrix(to[1], matrixA[1],matrixB);
 +      mul_fvector_fmatrix(to[2], matrixA[2],matrixB);
 +}
 +/* 4x4 matrix addition with 4x4 matrix */
 +DO_INLINE void add_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      VECADD(to[0], matrixA[0], matrixB[0]);
 +      VECADD(to[1], matrixA[1], matrixB[1]);
 +      VECADD(to[2], matrixA[2], matrixB[2]);
 +}
 +/* 4x4 matrix add-addition with 4x4 matrix */
 +DO_INLINE void addadd_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      VECADDADD(to[0], matrixA[0], matrixB[0]);
 +      VECADDADD(to[1], matrixA[1], matrixB[1]);
 +      VECADDADD(to[2], matrixA[2], matrixB[2]);
 +}
 +/* 4x4 matrix sub-addition with 4x4 matrix */
 +DO_INLINE void addsub_fmatrixS_fmatrixS(float to[3][3], float matrixA[3][3], float aS, float matrixB[3][3], float bS)
 +{
 +      VECADDSUBSS(to[0], matrixA[0], aS, matrixB[0], bS);
 +      VECADDSUBSS(to[1], matrixA[1], aS, matrixB[1], bS);
 +      VECADDSUBSS(to[2], matrixA[2], aS, matrixB[2], bS);
 +}
 +/* A -= B + C (4x4 matrix sub-addition with 4x4 matrix) */
 +DO_INLINE void subadd_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      VECSUBADD(to[0], matrixA[0], matrixB[0]);
 +      VECSUBADD(to[1], matrixA[1], matrixB[1]);
 +      VECSUBADD(to[2], matrixA[2], matrixB[2]);
 +}
 +/* A -= B*x + C*y (4x4 matrix sub-addition with 4x4 matrix) */
 +DO_INLINE void subadd_fmatrixS_fmatrixS(float to[3][3], float matrixA[3][3], float aS, float matrixB[3][3], float bS)
 +{
 +      VECSUBADDSS(to[0], matrixA[0], aS, matrixB[0], bS);
 +      VECSUBADDSS(to[1], matrixA[1], aS, matrixB[1], bS);
 +      VECSUBADDSS(to[2], matrixA[2], aS, matrixB[2], bS);
 +}
 +/* A = B - C (4x4 matrix subtraction with 4x4 matrix) */
 +DO_INLINE void sub_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      VECSUB(to[0], matrixA[0], matrixB[0]);
 +      VECSUB(to[1], matrixA[1], matrixB[1]);
 +      VECSUB(to[2], matrixA[2], matrixB[2]);
 +}
 +/* A += B - C (4x4 matrix add-subtraction with 4x4 matrix) */
 +DO_INLINE void addsub_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      VECADDSUB(to[0], matrixA[0], matrixB[0]);
 +      VECADDSUB(to[1], matrixA[1], matrixB[1]);
 +      VECADDSUB(to[2], matrixA[2], matrixB[2]);
 +}
 +/////////////////////////////////////////////////////////////////
 +// special functions
 +/////////////////////////////////////////////////////////////////
 +/* a vector multiplied and added to/by a 4x4 matrix */
 +DO_INLINE void muladd_fvector_fmatrix(float to[3], float from[3], float matrix[3][3])
 +{
 +      to[0] += matrix[0][0]*from[0] + matrix[1][0]*from[1] + matrix[2][0]*from[2];
 +      to[1] += matrix[0][1]*from[0] + matrix[1][1]*from[1] + matrix[2][1]*from[2];
 +      to[2] += matrix[0][2]*from[0] + matrix[1][2]*from[1] + matrix[2][2]*from[2];
 +}
 +/* 4x4 matrix multiplied and added  to/by a 4x4 matrix  and added to another 4x4 matrix */
 +DO_INLINE void muladd_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      muladd_fvector_fmatrix(to[0], matrixA[0],matrixB);
 +      muladd_fvector_fmatrix(to[1], matrixA[1],matrixB);
 +      muladd_fvector_fmatrix(to[2], matrixA[2],matrixB);
 +}
 +/* a vector multiplied and sub'd to/by a 4x4 matrix */
 +DO_INLINE void mulsub_fvector_fmatrix(float to[3], float from[3], float matrix[3][3])
 +{
 +      to[0] -= matrix[0][0]*from[0] + matrix[1][0]*from[1] + matrix[2][0]*from[2];
 +      to[1] -= matrix[0][1]*from[0] + matrix[1][1]*from[1] + matrix[2][1]*from[2];
 +      to[2] -= matrix[0][2]*from[0] + matrix[1][2]*from[1] + matrix[2][2]*from[2];
 +}
 +/* 4x4 matrix multiplied and sub'd  to/by a 4x4 matrix  and added to another 4x4 matrix */
 +DO_INLINE void mulsub_fmatrix_fmatrix(float to[3][3], float matrixA[3][3], float matrixB[3][3])
 +{
 +      mulsub_fvector_fmatrix(to[0], matrixA[0],matrixB);
 +      mulsub_fvector_fmatrix(to[1], matrixA[1],matrixB);
 +      mulsub_fvector_fmatrix(to[2], matrixA[2],matrixB);
 +}
 +/* 4x4 matrix multiplied+added by a vector */
 +/* STATUS: verified */
 +DO_INLINE void muladd_fmatrix_fvector(float to[3], float matrix[3][3], float from[3])
 +{
 +      to[0] += INPR(matrix[0],from);
 +      to[1] += INPR(matrix[1],from);
 +      to[2] += INPR(matrix[2],from);  
 +}
 +/* 4x4 matrix multiplied+sub'ed by a vector */
 +DO_INLINE void mulsub_fmatrix_fvector(float to[3], float matrix[3][3], float from[3])
 +{
 +      to[0] -= INPR(matrix[0],from);
 +      to[1] -= INPR(matrix[1],from);
 +      to[2] -= INPR(matrix[2],from);
 +}
 +/////////////////////////////////////////////////////////////////
 +
 +///////////////////////////
 +// SPARSE SYMMETRIC big matrix with 4x4 matrix entries
 +///////////////////////////
 +/* printf a big matrix on console: for debug output */
 +void print_bfmatrix(fmatrix3x3 *m3)
 +{
 +      unsigned int i = 0;
 +
 +      for(i = 0; i < m3[0].vcount + m3[0].scount; i++)
 +      {
 +              print_fmatrix(m3[i].m);
 +      }
 +}
 +/* create big matrix */
 +DO_INLINE fmatrix3x3 *create_bfmatrix(unsigned int verts, unsigned int springs)
 +{
 +      // TODO: check if memory allocation was successfull */
 +      fmatrix3x3 *temp = (fmatrix3x3 *)MEM_callocN (sizeof (fmatrix3x3) * (verts + springs), "cloth_implicit_alloc_matrix");
 +      temp[0].vcount = verts;
 +      temp[0].scount = springs;
 +      return temp;
 +}
 +/* delete big matrix */
 +DO_INLINE void del_bfmatrix(fmatrix3x3 *matrix)
 +{
 +      if (matrix != NULL)
 +      {
 +              MEM_freeN (matrix);
 +      }
 +}
 +/* copy big matrix */
 +DO_INLINE void cp_bfmatrix(fmatrix3x3 *to, fmatrix3x3 *from)
 +{     
 +      // TODO bounds checking 
 +      memcpy(to, from, sizeof(fmatrix3x3) * (from[0].vcount+from[0].scount) );
 +}
 +/* init the diagonal of big matrix */
 +// slow in parallel
 +DO_INLINE void initdiag_bfmatrix(fmatrix3x3 *matrix, float m3[3][3])
 +{
 +      unsigned int i,j;
 +      float tmatrix[3][3] = {{0,0,0},{0,0,0},{0,0,0}};
 +
 +      for(i = 0; i < matrix[0].vcount; i++)
 +      {               
 +              cp_fmatrix(matrix[i].m, m3); 
 +      }
 +      for(j = matrix[0].vcount; j < matrix[0].vcount+matrix[0].scount; j++)
 +      {
 +              cp_fmatrix(matrix[j].m, tmatrix); 
 +      }
 +}
 +/* init big matrix */
 +DO_INLINE void init_bfmatrix(fmatrix3x3 *matrix, float m3[3][3])
 +{
 +      unsigned int i;
 +
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              cp_fmatrix(matrix[i].m, m3); 
 +      }
 +}
 +/* multiply big matrix with scalar*/
 +DO_INLINE void mul_bfmatrix_S(fmatrix3x3 *matrix, float scalar)
 +{
 +      unsigned int i = 0;
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              mul_fmatrix_S(matrix[i].m, scalar);
 +      }
 +}
 +/* SPARSE SYMMETRIC multiply big matrix with long vector*/
 +/* STATUS: verified */
 +DO_INLINE void mul_bfmatrix_lfvector( float (*to)[3], fmatrix3x3 *from, lfVector *fLongVector)
 +{
 +      unsigned int i = 0;
 +      lfVector *temp = create_lfvector(from[0].vcount);
 +      
 +      zero_lfvector(to, from[0].vcount);
 +
 +#pragma omp parallel sections private(i)
 +      {
 +#pragma omp section
 +              {
 +                      for(i = from[0].vcount; i < from[0].vcount+from[0].scount; i++)
 +                      {
 +                              muladd_fmatrix_fvector(to[from[i].c], from[i].m, fLongVector[from[i].r]);
 +                      }
 +              }       
 +#pragma omp section
 +              {
 +                      for(i = 0; i < from[0].vcount+from[0].scount; i++)
 +                      {
 +                              muladd_fmatrix_fvector(temp[from[i].r], from[i].m, fLongVector[from[i].c]);
 +                      }
 +              }
 +      }
 +      add_lfvector_lfvector(to, to, temp, from[0].vcount);
 +      
 +      del_lfvector(temp);
 +      
 +      
 +}
 +/* SPARSE SYMMETRIC add big matrix with big matrix: A = B + C*/
 +DO_INLINE void add_bfmatrix_bfmatrix( fmatrix3x3 *to, fmatrix3x3 *from,  fmatrix3x3 *matrix)
 +{
 +      unsigned int i = 0;
 +
 +      /* process diagonal elements */
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              add_fmatrix_fmatrix(to[i].m, from[i].m, matrix[i].m);   
 +      }
 +
 +}
 +/* SPARSE SYMMETRIC add big matrix with big matrix: A += B + C */
 +DO_INLINE void addadd_bfmatrix_bfmatrix( fmatrix3x3 *to, fmatrix3x3 *from,  fmatrix3x3 *matrix)
 +{
 +      unsigned int i = 0;
 +
 +      /* process diagonal elements */
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              addadd_fmatrix_fmatrix(to[i].m, from[i].m, matrix[i].m);        
 +      }
 +
 +}
 +/* SPARSE SYMMETRIC subadd big matrix with big matrix: A -= B + C */
 +DO_INLINE void subadd_bfmatrix_bfmatrix( fmatrix3x3 *to, fmatrix3x3 *from,  fmatrix3x3 *matrix)
 +{
 +      unsigned int i = 0;
 +
 +      /* process diagonal elements */
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              subadd_fmatrix_fmatrix(to[i].m, from[i].m, matrix[i].m);        
 +      }
 +
 +}
 +/*  A = B - C (SPARSE SYMMETRIC sub big matrix with big matrix) */
 +DO_INLINE void sub_bfmatrix_bfmatrix( fmatrix3x3 *to, fmatrix3x3 *from,  fmatrix3x3 *matrix)
 +{
 +      unsigned int i = 0;
 +
 +      /* process diagonal elements */
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              sub_fmatrix_fmatrix(to[i].m, from[i].m, matrix[i].m);   
 +      }
 +
 +}
 +/* SPARSE SYMMETRIC sub big matrix with big matrix S (special constraint matrix with limited entries) */
 +DO_INLINE void sub_bfmatrix_Smatrix( fmatrix3x3 *to, fmatrix3x3 *from,  fmatrix3x3 *matrix)
 +{
 +      unsigned int i = 0;
 +
 +      /* process diagonal elements */
 +      for(i = 0; i < matrix[0].vcount; i++)
 +      {
 +              sub_fmatrix_fmatrix(to[matrix[i].c].m, from[matrix[i].c].m, matrix[i].m);       
 +      }
 +
 +}
 +/* A += B - C (SPARSE SYMMETRIC addsub big matrix with big matrix) */
 +DO_INLINE void addsub_bfmatrix_bfmatrix( fmatrix3x3 *to, fmatrix3x3 *from,  fmatrix3x3 *matrix)
 +{
 +      unsigned int i = 0;
 +
 +      /* process diagonal elements */
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              addsub_fmatrix_fmatrix(to[i].m, from[i].m, matrix[i].m);        
 +      }
 +
 +}
 +/* SPARSE SYMMETRIC sub big matrix with big matrix*/
 +/* A -= B * float + C * float --> for big matrix */
 +/* VERIFIED */
 +DO_INLINE void subadd_bfmatrixS_bfmatrixS( fmatrix3x3 *to, fmatrix3x3 *from, float aS,  fmatrix3x3 *matrix, float bS)
 +{
 +      unsigned int i = 0;
 +
 +      /* process diagonal elements */
 +      for(i = 0; i < matrix[0].vcount+matrix[0].scount; i++)
 +      {
 +              subadd_fmatrixS_fmatrixS(to[i].m, from[i].m, aS, matrix[i].m, bS);      
 +      }
 +
 +}
 +
 +///////////////////////////////////////////////////////////////////
 +// simulator start
 +///////////////////////////////////////////////////////////////////
 +static float I[3][3] = {{1,0,0},{0,1,0},{0,0,1}};
 +static float ZERO[3][3] = {{0,0,0}, {0,0,0}, {0,0,0}};
 +typedef struct Implicit_Data 
 +{
 +      lfVector *X, *V, *Xnew, *Vnew, *F, *B, *dV, *z;
 +      fmatrix3x3 *A, *dFdV, *dFdX, *S, *P, *Pinv, *bigI; 
 +} Implicit_Data;
 +
 +int implicit_init (Object *ob, ClothModifierData *clmd)
 +{
 +      unsigned int i = 0;
 +      unsigned int pinned = 0;
 +      Cloth *cloth = NULL;
 +      ClothVertex *verts = NULL;
 +      ClothSpring *spring = NULL;
 +      Implicit_Data *id = NULL;
 +      LinkNode *search = NULL;
 +
 +      // init memory guard
 +      // MEMORY_BASE.first = MEMORY_BASE.last = NULL;
 +
 +      cloth = (Cloth *)clmd->clothObject;
 +      verts = cloth->verts;
 +
 +      // create implicit base
 +      id = (Implicit_Data *)MEM_callocN (sizeof(Implicit_Data), "implicit vecmat");
 +      cloth->implicit = id;
 +
 +      /* process diagonal elements */         
 +      id->A = create_bfmatrix(cloth->numverts, cloth->numsprings);
 +      id->dFdV = create_bfmatrix(cloth->numverts, cloth->numsprings);
 +      id->dFdX = create_bfmatrix(cloth->numverts, cloth->numsprings);
 +      id->S = create_bfmatrix(cloth->numverts, 0);
 +      id->Pinv = create_bfmatrix(cloth->numverts, cloth->numsprings);
 +      id->P = create_bfmatrix(cloth->numverts, cloth->numsprings);
 +      id->bigI = create_bfmatrix(cloth->numverts, cloth->numsprings); // TODO 0 springs
 +      id->X = create_lfvector(cloth->numverts);
 +      id->Xnew = create_lfvector(cloth->numverts);
 +      id->V = create_lfvector(cloth->numverts);
 +      id->Vnew = create_lfvector(cloth->numverts);
 +      id->F = create_lfvector(cloth->numverts);
 +      id->B = create_lfvector(cloth->numverts);
 +      id->dV = create_lfvector(cloth->numverts);
 +      id->z = create_lfvector(cloth->numverts);
 +      
 +      for(i=0;i<cloth->numverts;i++) 
 +      {
 +              id->A[i].r = id->A[i].c = id->dFdV[i].r = id->dFdV[i].c = id->dFdX[i].r = id->dFdX[i].c = id->P[i].c = id->P[i].r = id->Pinv[i].c = id->Pinv[i].r = id->bigI[i].c = id->bigI[i].r = i;
 +
 +              if(verts [i].goal >= SOFTGOALSNAP)
 +              {
 +                      id->S[pinned].pinned = 1;
 +                      id->S[pinned].c = id->S[pinned].r = i;
 +                      pinned++;
 +              }
 +      }
 +
 +      // S is special and needs specific vcount and scount
 +      id->S[0].vcount = pinned; id->S[0].scount = 0;
 +
 +      // init springs */
 +      search = cloth->springs;
 +      for(i=0;i<cloth->numsprings;i++) 
 +      {
 +              spring = search->link;
 +              
 +              // dFdV_start[i].r = big_I[i].r = big_zero[i].r = 
 +              id->A[i+cloth->numverts].r = id->dFdV[i+cloth->numverts].r = id->dFdX[i+cloth->numverts].r = 
 +                              id->P[i+cloth->numverts].r = id->Pinv[i+cloth->numverts].r = id->bigI[i+cloth->numverts].r = spring->ij;
 +
 +              // dFdV_start[i].c = big_I[i].c = big_zero[i].c = 
 +              id->A[i+cloth->numverts].c = id->dFdV[i+cloth->numverts].c = id->dFdX[i+cloth->numverts].c = 
 +                              id->P[i+cloth->numverts].c = id->Pinv[i+cloth->numverts].c = id->bigI[i+cloth->numverts].c = spring->kl;
 +
 +              spring->matrix_index = i + cloth->numverts;
 +              
 +              search = search->next;
 +      }
 +
 +      for(i = 0; i < cloth->numverts; i++)
 +      {               
 +              VECCOPY(id->X[i], cloth->x[i]);
 +      }
 +
 +      return 1;
 +}
 +int   implicit_free (ClothModifierData *clmd)
 +{
 +      Implicit_Data *id;
 +      Cloth *cloth;
 +      cloth = (Cloth *)clmd->clothObject;
 +
 +      if(cloth)
 +      {
 +              id = cloth->implicit;
 +
 +              if(id)
 +              {
 +                      del_bfmatrix(id->A);
 +                      del_bfmatrix(id->dFdV);
 +                      del_bfmatrix(id->dFdX);
 +                      del_bfmatrix(id->S);
 +                      del_bfmatrix(id->P);
 +                      del_bfmatrix(id->Pinv);
 +                      del_bfmatrix(id->bigI);
 +
 +                      del_lfvector(id->X);
 +                      del_lfvector(id->Xnew);
 +                      del_lfvector(id->V);
 +                      del_lfvector(id->Vnew);
 +                      del_lfvector(id->F);
 +                      del_lfvector(id->B);
 +                      del_lfvector(id->dV);
 +                      del_lfvector(id->z);
 +
 +                      MEM_freeN(id);
 +              }
 +      }
 +
 +      return 1;
 +}
 +
 +void cloth_bending_mode(ClothModifierData *clmd, int enabled)
 +{
 +      Cloth *cloth = clmd->clothObject;
 +      Implicit_Data *id;
 +      
 +      if(cloth)
 +      {
 +              id = cloth->implicit;
 +              
 +              if(id)
 +              {
 +                      if(enabled)
 +                      {
 +                              cloth->numsprings = cloth->numspringssave;
 +                      }
 +                      else
 +                      {
 +                              cloth->numsprings = cloth->numothersprings;
 +                      }
 +                      
 +                      id->A[0].scount = id->dFdV[0].scount = id->dFdX[0].scount = id->P[0].scount = id->Pinv[0].scount = id->bigI[0].scount = cloth->numsprings;
 +              }       
 +      }
 +}
 +
 +DO_INLINE float fb(float length, float L)
 +{
 +      float x = length/L;
 +      return (-11.541f*pow(x,4)+34.193f*pow(x,3)-39.083f*pow(x,2)+23.116f*x-9.713f);
 +}
 +
 +DO_INLINE float fbderiv(float length, float L)
 +{
 +      float x = length/L;
 +
 +      return (-46.164f*pow(x,3)+102.579f*pow(x,2)-78.166f*x+23.116f);
 +}
 +
 +DO_INLINE float fbstar(float length, float L, float kb, float cb)
 +{
 +      float tempfb = kb * fb(length, L);
 +
 +      float fbstar = cb * (length - L);
 +
 +      if(tempfb < fbstar)
 +              return fbstar;
 +      else
 +              return tempfb;          
 +}
 +
 +DO_INLINE float fbstar_jacobi(float length, float L, float kb, float cb)
 +{
 +      float tempfb = kb * fb(length, L);
 +      float fbstar = cb * (length - L);
 +
 +      if(tempfb < fbstar)
 +      {               
 +              return cb;
 +      }
 +      else
 +      {
 +              return kb * fbderiv(length, L); 
 +      }       
 +}
 +
 +DO_INLINE void filter(lfVector *V, fmatrix3x3 *S)
 +{
 +      unsigned int i=0;
 +
 +      for(i=0;i<S[0].vcount;i++)
 +      {
 +              mul_fvector_fmatrix(V[S[i].r], V[S[i].r], S[i].m);
 +      }
 +}
 +
 +int cg_filtered(lfVector *dv, fmatrix3x3 *lA, lfVector *lB, lfVector *z, fmatrix3x3 *S)
 +{
 +      // Solves for unknown X in equation AX=B
 +      unsigned int conjgrad_loopcount=0, conjgrad_looplimit=100;
 +      float conjgrad_epsilon=0.0001f, conjgrad_lasterror=0;
 +      lfVector *q, *d, *tmp, *r; 
 +      float s, starget, a, s_prev;
 +      unsigned int numverts = lA[0].vcount;
 +      q = create_lfvector(numverts);
 +      d = create_lfvector(numverts);
 +      tmp = create_lfvector(numverts);
 +      r = create_lfvector(numverts);
 +
 +      // zero_lfvector(dv, CLOTHPARTICLES);
 +      filter(dv, S);
 +
 +      add_lfvector_lfvector(dv, dv, z, numverts);
 +
 +      // r = B - Mul(tmp,A,X);    // just use B if X known to be zero
 +      cp_lfvector(r, lB, numverts);
 +      mul_bfmatrix_lfvector(tmp, lA, dv);
 +      sub_lfvector_lfvector(r, r, tmp, numverts);
 +
 +      filter(r,S);
 +
 +      cp_lfvector(d, r, numverts);
 +
 +      s = dot_lfvector(r, r, numverts);
 +      starget = s * sqrt(conjgrad_epsilon);
 +
 +      while((s>starget && conjgrad_loopcount < conjgrad_looplimit))
 +      {       
 +              // Mul(q,A,d); // q = A*d;
 +              mul_bfmatrix_lfvector(q, lA, d);
 +
 +              filter(q,S);
 +
 +              a = s/dot_lfvector(d, q, numverts);
 +
 +              // X = X + d*a;
 +              add_lfvector_lfvectorS(dv, dv, d, a, numverts);
 +
 +              // r = r - q*a;
 +              sub_lfvector_lfvectorS(r, r, q, a, numverts);
 +
 +              s_prev = s;
 +              s = dot_lfvector(r, r, numverts);
 +
 +              //d = r+d*(s/s_prev);
 +              add_lfvector_lfvectorS(d, r, d, (s/s_prev), numverts);
 +
 +              filter(d,S);
 +
 +              conjgrad_loopcount++;
 +      }
 +      // itend();
 +      // printf("cg_filtered time: %f\n", (float)itval());
 +      
 +      conjgrad_lasterror = s;
 +
 +      del_lfvector(q);
 +      del_lfvector(d);
 +      del_lfvector(tmp);
 +      del_lfvector(r);
 +      
 +      // printf("W/O conjgrad_loopcount: %d\n", conjgrad_loopcount);
 +
 +      return conjgrad_loopcount<conjgrad_looplimit;  // true means we reached desired accuracy in given time - ie stable
 +}
 +
 +// block diagonalizer
 +DO_INLINE void BuildPPinv(fmatrix3x3 *lA, fmatrix3x3 *P, fmatrix3x3 *Pinv)
 +{
 +      unsigned int i = 0;
 +      
 +      // Take only the diagonal blocks of A
 +// #pragma omp parallel for private(i)
 +      for(i = 0; i<lA[0].vcount; i++)
 +      {
 +              // block diagonalizer
 +              cp_fmatrix(P[i].m, lA[i].m);
 +              inverse_fmatrix(Pinv[i].m, P[i].m);
 +              
 +      }
 +}
 +
 +// version 1.3
 +int cg_filtered_pre(lfVector *dv, fmatrix3x3 *lA, lfVector *lB, lfVector *z, fmatrix3x3 *S, fmatrix3x3 *P, fmatrix3x3 *Pinv)
 +{
 +      unsigned int numverts = lA[0].vcount, iterations = 0, conjgrad_looplimit=100;
 +      float delta0 = 0, deltaNew = 0, deltaOld = 0, alpha = 0;
 +      float conjgrad_epsilon=0.0001; // 0.2 is dt for steps=5
 +      lfVector *r = create_lfvector(numverts);
 +      lfVector *p = create_lfvector(numverts);
 +      lfVector *s = create_lfvector(numverts);
 +      lfVector *h = create_lfvector(numverts);
 +      
 +      BuildPPinv(lA, P, Pinv);
 +      
 +      filter(dv, S);
 +      add_lfvector_lfvector(dv, dv, z, numverts);
 +      
 +      mul_bfmatrix_lfvector(r, lA, dv);
 +      sub_lfvector_lfvector(r, lB, r, numverts);
 +      filter(r, S);
 +      
 +      mul_bfmatrix_lfvector(p, Pinv, r);
 +      filter(p, S);
 +      
 +      deltaNew = dot_lfvector(r, p, numverts);
 +      
 +      delta0 = deltaNew * sqrt(conjgrad_epsilon);
 +      
 +      // itstart();
 +      
 +      while ((deltaNew > delta0) && (iterations < conjgrad_looplimit))
 +      {
 +              iterations++;
 +              
 +              mul_bfmatrix_lfvector(s, lA, p);
 +              filter(s, S);
 +              
 +              alpha = deltaNew / dot_lfvector(p, s, numverts);
 +              
 +              add_lfvector_lfvectorS(dv, dv, p, alpha, numverts);
 +              
 +              add_lfvector_lfvectorS(r, r, s, -alpha, numverts);
 +              
 +              mul_bfmatrix_lfvector(h, Pinv, r);
 +              filter(h, S);
 +              
 +              deltaOld = deltaNew;
 +              
 +              deltaNew = dot_lfvector(r, h, numverts);
 +              
 +              add_lfvector_lfvectorS(p, h, p, deltaNew / deltaOld, numverts);
 +              filter(p, S);
 +              
 +      }
 +      
 +      // itend();
 +      // printf("cg_filtered_pre time: %f\n", (float)itval());
 +      
 +      del_lfvector(h);
 +      del_lfvector(s);
 +      del_lfvector(p);
 +      del_lfvector(r);
 +      
 +      // printf("iterations: %d\n", iterations);
 +      
 +      return iterations<conjgrad_looplimit;
 +}
 +
 +// outer product is NOT cross product!!!
 +DO_INLINE void dfdx_spring_type1(float to[3][3], float dir[3],float length,float L,float k)
 +{
 +      // dir is unit length direction, rest is spring's restlength, k is spring constant.
 +      // return  (outerprod(dir,dir)*k + (I - outerprod(dir,dir))*(k - ((k*L)/length)));
 +      float temp[3][3];
 +      mul_fvectorT_fvector(temp, dir, dir);
 +      sub_fmatrix_fmatrix(to, I, temp);
 +      mul_fmatrix_S(to, k* (1.0f-(L/length)));
 +      mul_fmatrix_S(temp, k);
 +      add_fmatrix_fmatrix(to, temp, to);
 +}
 +
 +DO_INLINE void dfdx_spring_type2(float to[3][3], float dir[3],float length,float L,float k, float cb)
 +{
 +      // return  outerprod(dir,dir)*fbstar_jacobi(length, L, k, cb);
 +      mul_fvectorT_fvectorS(to, dir, dir, fbstar_jacobi(length, L, k, cb));
 +}
 +
 +DO_INLINE void dfdv_damp(float to[3][3], float dir[3], float damping)
 +{
 +      // derivative of force wrt velocity.  
 +      // return outerprod(dir,dir) * damping;
 +      mul_fvectorT_fvectorS(to, dir, dir, damping);
 +}
 +
 +DO_INLINE void dfdx_spring(float to[3][3],  float dir[3],float length,float L,float k)
 +{
 +      // dir is unit length direction, rest is spring's restlength, k is spring constant.
 +      //return  ( (I-outerprod(dir,dir))*Min(1.0f,rest/length) - I) * -k;
 +      mul_fvectorT_fvector(to, dir, dir);
 +      sub_fmatrix_fmatrix(to, I, to);
 +      mul_fmatrix_S(to, (((L/length)> 1.0f) ? (1.0f): (L/length))); 
 +      sub_fmatrix_fmatrix(to, to, I);
 +      mul_fmatrix_S(to, -k);
 +}
 +
 +DO_INLINE void dfdx_damp(float to[3][3],  float dir[3],float length,const float vel[3],float rest,float damping)
 +{
 +      // inner spring damping   vel is the relative velocity  of the endpoints.  
 +      //      return (I-outerprod(dir,dir)) * (-damping * -(dot(dir,vel)/Max(length,rest)));
 +      mul_fvectorT_fvector(to, dir, dir);
 +      sub_fmatrix_fmatrix(to, I, to);
 +      mul_fmatrix_S(to,  (-damping * -(INPR(dir,vel)/MAX2(length,rest)))); 
 +
 +}
 +
 +DO_INLINE void cloth_calc_spring_force(ClothModifierData *clmd, ClothSpring *s, lfVector *lF, lfVector *X, lfVector *V, fmatrix3x3 *dFdV, fmatrix3x3 *dFdX)
 +{
 +      float extent[3];
 +      float length = 0;
 +      float dir[3] = {0,0,0};
 +      float vel[3];
 +      float k = 0.0f;
 +      float L = s->restlen;
 +      float cb = clmd->sim_parms->structural;
 +
 +      float nullf[3] = {0,0,0};
 +      float stretch_force[3] = {0,0,0};
 +      float bending_force[3] = {0,0,0};
 +      float damping_force[3] = {0,0,0};
 +      float nulldfdx[3][3]={ {0,0,0}, {0,0,0}, {0,0,0}};
 +      
 +      VECCOPY(s->f, nullf);
 +      cp_fmatrix(s->dfdx, nulldfdx);
 +      cp_fmatrix(s->dfdv, nulldfdx);
 +
 +      // calculate elonglation
 +      VECSUB(extent, X[s->kl], X[s->ij]);
 +      VECSUB(vel, V[s->kl], V[s->ij]);
 +      length = sqrt(INPR(extent, extent));
 +      
 +      s->flags &= ~CLOTH_SPRING_FLAG_NEEDED;
 +      
 +      if(length > ABS(ALMOST_ZERO))
 +      {
 +              /*
 +              if(length>L)
 +              {
 +              if((clmd->sim_parms->flags & CSIMSETT_FLAG_TEARING_ENABLED) 
 +              && ((((length-L)*100.0f/L) > clmd->sim_parms->maxspringlen))) // cut spring!
 +              {
 +              s->flags |= CSPRING_FLAG_DEACTIVATE;
 +              return;
 +      }
 +      } 
 +              */
 +              mul_fvector_S(dir, extent, 1.0f/length);
 +      }
 +      else    
 +      {
 +              mul_fvector_S(dir, extent, 0.0f);
 +      }
 +      
 +      
 +      // calculate force of structural + shear springs
 +      if(s->type != CLOTH_SPRING_TYPE_BENDING)
 +      {
 +              if(length > L) // only on elonglation
 +              {
 +                      s->flags |= CLOTH_SPRING_FLAG_NEEDED;
 +
 +                      k = clmd->sim_parms->structural;        
 +
 +                      mul_fvector_S(stretch_force, dir, (k*(length-L))); 
 +
 +                      VECADD(s->f, s->f, stretch_force);
 +
 +                      // Ascher & Boxman, p.21: Damping only during elonglation
 +                      mul_fvector_S(damping_force, extent, clmd->sim_parms->Cdis * ((INPR(vel,extent)/length))); 
 +                      VECADD(s->f, s->f, damping_force);
 +
 +                      dfdx_spring_type1(s->dfdx, dir,length,L,k);
 +
 +                      dfdv_damp(s->dfdv, dir,clmd->sim_parms->Cdis);
 +              }
 +      }
 +      else // calculate force of bending springs
 +      {
 +              if(length < L)
 +              {
 +                      // clmd->sim_parms->flags |= CLOTH_SIMSETTINGS_FLAG_BIG_FORCE;
 +                      
 +                      s->flags |= CLOTH_SPRING_FLAG_NEEDED;
 +                      
 +                      k = clmd->sim_parms->bending;   
 +
 +                      mul_fvector_S(bending_force, dir, fbstar(length, L, k, cb));
 +                      VECADD(s->f, s->f, bending_force);
 +                      
 +                      if(INPR(bending_force,bending_force) > 0.13*0.13)
 +                      {
 +                              clmd->sim_parms->flags |= CLOTH_SIMSETTINGS_FLAG_BIG_FORCE;
 +                      }
 +                      
 +                      dfdx_spring_type2(s->dfdx, dir,length,L,k, cb);
 +              }
 +      }
 +}
 +
 +DO_INLINE int cloth_apply_spring_force(ClothModifierData *clmd, ClothSpring *s, lfVector *lF, lfVector *X, lfVector *V, fmatrix3x3 *dFdV, fmatrix3x3 *dFdX)
 +{
 +      if(s->flags & CLOTH_SPRING_FLAG_NEEDED)
 +      {
 +              VECADD(lF[s->ij], lF[s->ij], s->f);
 +              VECSUB(lF[s->kl], lF[s->kl], s->f);     
 +                      
 +              if(s->type != CLOTH_SPRING_TYPE_BENDING)
 +              {
 +                      sub_fmatrix_fmatrix(dFdV[s->ij].m, dFdV[s->ij].m, s->dfdv);
 +                      sub_fmatrix_fmatrix(dFdV[s->kl].m, dFdV[s->kl].m, s->dfdv);
 +                      add_fmatrix_fmatrix(dFdV[s->matrix_index].m, dFdV[s->matrix_index].m, s->dfdv); 
 +              }
 +              else if(!(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_BIG_FORCE))
 +                      return 0;
 +              
 +              sub_fmatrix_fmatrix(dFdX[s->ij].m, dFdX[s->ij].m, s->dfdx);
 +              sub_fmatrix_fmatrix(dFdX[s->kl].m, dFdX[s->kl].m, s->dfdx);
 +
 +              add_fmatrix_fmatrix(dFdX[s->matrix_index].m, dFdX[s->matrix_index].m, s->dfdx);
 +      }
 +      
 +      return 1;
 +}
 +
 +DO_INLINE void calculateTriangleNormal(float to[3], lfVector *X, MFace mface)
 +{
 +      float v1[3], v2[3];
 +
 +      VECSUB(v1, X[mface.v2], X[mface.v1]);
 +      VECSUB(v2, X[mface.v3], X[mface.v1]);
 +      cross_fvector(to, v1, v2);
 +}
 +
 +DO_INLINE void calculatQuadNormal(float to[3], lfVector *X, MFace mface)
 +{
 +      float temp = CalcNormFloat4(X[mface.v1],X[mface.v2],X[mface.v3],X[mface.v4],to);
 +      mul_fvector_S(to, to, temp);
 +}
 +
 +void calculateWeightedVertexNormal(ClothModifierData *clmd, MFace *mfaces, float to[3], int index, lfVector *X)
 +{
 +      float temp[3]; 
 +      int i;
 +      Cloth *cloth = clmd->clothObject;
 +
 +      for(i = 0; i < cloth->numfaces; i++)
 +      {
 +              // check if this triangle contains the selected vertex
 +              if(mfaces[i].v1 == index || mfaces[i].v2 == index || mfaces[i].v3 == index || mfaces[i].v4 == index)
 +              {
 +                      calculatQuadNormal(temp, X, mfaces[i]);
 +                      VECADD(to, to, temp);
 +              }
 +      }
 +}
 +float calculateVertexWindForce(float wind[3], float vertexnormal[3])  
 +{
 +      return fabs(INPR(wind, vertexnormal) * 0.5f);
 +}
 +
 +DO_INLINE void calc_triangle_force(ClothModifierData *clmd, MFace mface, lfVector *F, lfVector *X, lfVector *V, fmatrix3x3 *dFdV, fmatrix3x3 *dFdX, ListBase *effectors)
 +{     
 +
 +}
 +
 +void cloth_calc_force(ClothModifierData *clmd, lfVector *lF, lfVector *lX, lfVector *lV, fmatrix3x3 *dFdV, fmatrix3x3 *dFdX, ListBase *effectors, float time)
 +{
 +      /* Collect forces and derivatives:  F,dFdX,dFdV */
 +      Cloth           *cloth          = clmd->clothObject;
 +      unsigned int    i               = 0;
 +      float           spring_air      = clmd->sim_parms->Cvi * 0.01f; /* viscosity of air scaled in percent */
 +      float           gravity[3];
 +      float           tm2[3][3]       = {{-spring_air,0,0}, {0,-spring_air,0},{0,0,-spring_air}};
 +      ClothVertex *verts = cloth->verts;
 +      MFace           *mfaces         = cloth->mfaces;
 +      float wind_normalized[3];
 +      unsigned int numverts = cloth->numverts;
 +      float auxvect[3], velgoal[3], tvect[3];
 +      float kd, ks;
 +      LinkNode *search = cloth->springs;
 +
 +      VECCOPY(gravity, clmd->sim_parms->gravity);
 +      mul_fvector_S(gravity, gravity, 0.001f); /* scale gravity force */
 +
 +      /* set dFdX jacobi matrix to zero */
 +      init_bfmatrix(dFdX, ZERO);
 +      /* set dFdX jacobi matrix diagonal entries to -spring_air */ 
 +      initdiag_bfmatrix(dFdV, tm2);
 +
 +      init_lfvector(lF, gravity, numverts);
 +
 +      submul_lfvectorS(lF, lV, spring_air, numverts);
 +
 +      /* do goal stuff */
 +      if(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) 
 +      {       
 +              for(i = 0; i < numverts; i++)
 +              {                       
 +                      if(verts [i].goal < SOFTGOALSNAP)
 +                      {                       
 +                              // current_position = xold + t * (newposition - xold)
 +                              VECSUB(tvect, cloth->xconst[i], cloth->xold[i]);
 +                              mul_fvector_S(tvect, tvect, time);
 +                              VECADD(tvect, tvect, cloth->xold[i]);
 +
 +                              VECSUB(auxvect, tvect, lX[i]);
 +                              ks  = 1.0f/(1.0f- verts [i].goal*clmd->sim_parms->goalspring)-1.0f ;
 +                              VECADDS(lF[i], lF[i], auxvect, -ks);
 +
 +                              // calulate damping forces generated by goals                           
 +                              VECSUB(velgoal, cloth->xold[i], cloth->xconst[i]);
 +                              kd =  clmd->sim_parms->goalfrict * 0.01f; // friction force scale taken from SB
 +                              VECSUBADDSS(lF[i], velgoal, kd, lV[i], kd);
 +                              
 +                      }
 +              }       
 +      }
 +
 +      /* handle external forces like wind */
 +      if(effectors)
 +      {
 +              float speed[3] = {0.0f, 0.0f,0.0f};
 +              float force[3]= {0.0f, 0.0f, 0.0f};
 +              
 +#pragma omp parallel for private (i) shared(lF)
 +              for(i = 0; i < cloth->numverts; i++)
 +              {
 +                      float vertexnormal[3]={0,0,0};
 +                      float fieldfactor = 1000.0f, windfactor  = 250.0f; // from sb
 +                      
 +                      pdDoEffectors(effectors, lX[i], force, speed, (float)G.scene->r.cfra, 0.0f, PE_WIND_AS_SPEED);          
 +                      
 +                      // TODO apply forcefields here
 +                      VECADDS(lF[i], lF[i], force, fieldfactor*0.01f);
 +
 +                      VECCOPY(wind_normalized, speed);
 +                      Normalize(wind_normalized);
 +                      
 +                      calculateWeightedVertexNormal(clmd, mfaces, vertexnormal, i, lX);
 +                      VECADDS(lF[i], lF[i], wind_normalized, -calculateVertexWindForce(speed, vertexnormal));
 +              }
 +      }
 +      
 +      // calculate spring forces
 +      search = cloth->springs;
 +      while(search)
 +      {
 +              // only handle active springs
 +              // if(((clmd->sim_parms->flags & CSIMSETT_FLAG_TEARING_ENABLED) && !(springs[i].flags & CSPRING_FLAG_DEACTIVATE))|| !(clmd->sim_parms->flags & CSIMSETT_FLAG_TEARING_ENABLED)){}
 +              cloth_calc_spring_force(clmd, search->link, lF, lX, lV, dFdV, dFdX);
 +
 +              search = search->next;
 +      }
 +      
 +      if(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_BIG_FORCE)
 +      {       
 +              if(cloth->numspringssave != cloth->numsprings)
 +              {
 +                      cloth_bending_mode(clmd, 1);
 +              }
 +      }
 +      else
 +      {
 +              if(cloth->numspringssave == cloth->numsprings)
 +              {
 +                      cloth_bending_mode(clmd, 0);
 +              }
 +      }
 +      
 +      // apply spring forces
 +      search = cloth->springs;
 +      while(search)
 +      {
 +              // only handle active springs
 +              // if(((clmd->sim_parms->flags & CSIMSETT_FLAG_TEARING_ENABLED) && !(springs[i].flags & CSPRING_FLAG_DEACTIVATE))|| !(clmd->sim_parms->flags & CSIMSETT_FLAG_TEARING_ENABLED))  
 +              if(!cloth_apply_spring_force(clmd, search->link, lF, lX, lV, dFdV, dFdX))
 +                      break;
 +              search = search->next;
 +      }
 +      
 +      clmd->sim_parms->flags &= ~CLOTH_SIMSETTINGS_FLAG_BIG_FORCE;
 +      
 +}
 +
 +void simulate_implicit_euler(lfVector *Vnew, lfVector *lX, lfVector *lV, lfVector *lF, fmatrix3x3 *dFdV, fmatrix3x3 *dFdX, float dt, fmatrix3x3 *A, lfVector *B, lfVector *dV, fmatrix3x3 *S, lfVector *z, fmatrix3x3 *P, fmatrix3x3 *Pinv)
 +{
 +      unsigned int numverts = dFdV[0].vcount;
 +
 +      lfVector *dFdXmV = create_lfvector(numverts);
 +      initdiag_bfmatrix(A, I);
 +      zero_lfvector(dV, numverts);
 +
 +      subadd_bfmatrixS_bfmatrixS(A, dFdV, dt, dFdX, (dt*dt));
 +
 +      mul_bfmatrix_lfvector(dFdXmV, dFdX, lV);
 +
 +      add_lfvectorS_lfvectorS(B, lF, dt, dFdXmV, (dt*dt), numverts);
 +      
 +      // itstart();
 +      
 +      cg_filtered(dV, A, B, z, S); /* conjugate gradient algorithm to solve Ax=b */
 +      
 +      // TODO: if anyone finds a way to correct this function =>
 +      // explodes with stiffness = 3000 and 16k verts + pinned at 2 corners
 +      // cg_filtered_pre(dV, A, B, z, S, P, Pinv);
 +      
 +      // itend();
 +      // printf("cg_filtered calc time: %f\n", (float)itval());
 +      
 +      // advance velocities
 +      add_lfvector_lfvector(Vnew, lV, dV, numverts);
 +      
 +
 +      del_lfvector(dFdXmV);
 +}
 +
 +int implicit_solver (Object *ob, float frame, ClothModifierData *clmd, ListBase *effectors)
 +{             
 +      unsigned int i=0;
 +      float step=0.0f, tf=1.0f;
 +      Cloth *cloth = clmd->clothObject;
 +      ClothVertex *verts = cloth->verts;
 +      unsigned int numverts = cloth->numverts;
 +      float dt = 1.0f / clmd->sim_parms->stepsPerFrame;
 +      Implicit_Data *id = cloth->implicit;
 +      int result = 0;
 +      float force = 0, lastforce = 0;
 +      // lfVector *dx;
 +      
 +      if(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) /* do goal stuff */
 +      {
 +              for(i = 0; i < numverts; i++)
 +              {                       
 +                      // update velocities with constrained velocities from pinned verts
 +                      if(verts [i].goal >= SOFTGOALSNAP)
 +                      {                       
 +                              VECSUB(id->V[i], cloth->xconst[i], cloth->xold[i]);
 +                              // VecMulf(id->V[i], 1.0 / dt);
 +                      }
 +              }       
 +      }
 +
 +      while(step < tf)
 +      {               
 +              effectors= pdInitEffectors(ob,NULL);
 +              
 +              // calculate 
 +              cloth_calc_force(clmd, id->F, id->X, id->V, id->dFdV, id->dFdX, effectors, step );
 +              
 +              // check for sleeping
 +              // if(!(clmd->coll_parms->flags & CLOTH_SIMSETTINGS_FLAG_SLEEP))
 +              {
 +                      simulate_implicit_euler(id->Vnew, id->X, id->V, id->F, id->dFdV, id->dFdX, dt, id->A, id->B, id->dV, id->S, id->z, id->P, id->Pinv);
 +              
 +                      add_lfvector_lfvectorS(id->Xnew, id->X, id->Vnew, dt, numverts);
 +              }
 +              /*
 +              dx = create_lfvector(numverts);
 +              sub_lfvector_lfvector(dx, id->Xnew, id->X, numverts);
 +              force = dot_lfvector(dx, dx, numverts);
 +              del_lfvector(dx);
 +              
 +              if((force < 0.00001) && (lastforce >= force))
 +              clmd->coll_parms->flags |= CLOTH_SIMSETTINGS_FLAG_SLEEP;
 +              else if((lastforce*2 < force))
 +              clmd->coll_parms->flags &= ~CLOTH_SIMSETTINGS_FLAG_SLEEP;
 +              */
 +              lastforce = force;
 +              
 +              if(clmd->coll_parms->flags & CLOTH_COLLISIONSETTINGS_FLAG_ENABLED)
 +              {
 +                      // collisions 
 +                      // itstart();
 +                      
 +                      // update verts to current positions
 +                      for(i = 0; i < numverts; i++)
 +                      {               
 +                              if(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) /* do goal stuff */
 +                              {                       
 +                                      if(verts [i].goal >= SOFTGOALSNAP)
 +                                      {                       
 +                                              float tvect[3] = {.0,.0,.0};
 +                                              // VECSUB(tvect, id->Xnew[i], verts[i].xold);
 +                                              mul_fvector_S(tvect, id->V[i], step+dt);
 +                                              VECADD(tvect, tvect, cloth->xold[i]);
 +                                              VECCOPY(id->Xnew[i], tvect);
 +                                      }
 +                                              
 +                              }
 +                              
 +                              VECCOPY(cloth->current_x[i], id->Xnew[i]);
 +                              VECSUB(cloth->current_v[i], cloth->current_x[i], cloth->current_xold[i]);
 +                              VECCOPY(cloth->v[i], cloth->current_v[i]);
 +                      }
 +                      
 +                      // call collision function
 +                      result = cloth_bvh_objcollision(clmd, step + dt, step, dt);
 +                      
 +                      // copy corrected positions back to simulation
 +                      if(result)
 +                      {
 +                              memcpy(cloth->current_xold, cloth->current_x, sizeof(lfVector) * numverts);
 +                              memcpy(id->Xnew, cloth->current_x, sizeof(lfVector) * numverts);
 +                              
 +                              for(i = 0; i < numverts; i++)
 +                              {       
 +                                      VECCOPY(id->Vnew[i], cloth->current_v[i]);
 +                                      VecMulf(id->Vnew[i], 1.0f / dt);
 +                              }
 +                      }
 +                      else
 +                      {
 +                              memcpy(cloth->current_xold, id->Xnew, sizeof(lfVector) * numverts);
 +                      }
 +                      
 +                      // X = Xnew;
 +                      cp_lfvector(id->X, id->Xnew, numverts);
 +                      
 +                      // if there were collisions, advance the velocity from v_n+1/2 to v_n+1
 +                      if(result)
 +                      {
 +                              // V = Vnew;
 +                              cp_lfvector(id->V, id->Vnew, numverts);
 +                              
 +                              // calculate 
 +                              cloth_calc_force(clmd, id->F, id->X, id->V, id->dFdV, id->dFdX, effectors, step);       
 +                              simulate_implicit_euler(id->Vnew, id->X, id->V, id->F, id->dFdV, id->dFdX, dt / 2.0f, id->A, id->B, id->dV, id->S, id->z, id->P, id->Pinv);
 +                      }
 +                      
 +              }
 +              else
 +              {
 +                      // X = Xnew;
 +                      cp_lfvector(id->X, id->Xnew, numverts);
 +              }
 +              
 +              // itend();
 +              // printf("collision time: %f\n", (float)itval());
 +              
 +              // V = Vnew;
 +              cp_lfvector(id->V, id->Vnew, numverts);
 +              
 +              step += dt;
 +
 +              if(effectors) pdEndEffectors(effectors);
 +      }
 +      
 +      if(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL)
 +      {
 +              for(i = 0; i < numverts; i++)
 +              {
 +                      if(verts [i].goal < SOFTGOALSNAP)
 +                      {
 +                              VECCOPY(cloth->current_xold[i], id->X[i]);
 +                              VECCOPY(cloth->x[i], id->X[i]);
 +                      }
 +                      else
 +                      {
 +                              VECCOPY(cloth->current_xold[i], cloth->xconst[i]);
 +                              VECCOPY(cloth->x[i], cloth->xconst[i]);
 +                      }
 +              }
 +      }
 +      else
 +      {
 +              memcpy(cloth->current_xold, id->X, sizeof(lfVector) * numverts);
 +              memcpy(cloth->x, id->X, sizeof(lfVector) * numverts);
 +      }
 +      
 +      memcpy(cloth->v, id->V, sizeof(lfVector) * numverts);
 +      
 +      return 1;
 +}
 +
 +void implicit_set_positions (ClothModifierData *clmd)
 +{             
 +      Cloth *cloth = clmd->clothObject;
 +      unsigned int numverts = cloth->numverts;
 +      Implicit_Data *id = cloth->implicit;
 +      
 +      memcpy(id->X, cloth->x, sizeof(lfVector) * numverts);   
 +      memcpy(id->V, cloth->v, sizeof(lfVector) * numverts);   
 +}
 +
 +
 +int collisions_collision_response_static(ClothModifierData *clmd, ClothModifierData *coll_clmd)
 +{
 +      /*
 +      unsigned int i = 0;
 +      int result = 0;
 +      LinkNode *search = NULL;
 +      CollPair *collpair = NULL;
 +      Cloth *cloth1, *cloth2;
 +      float w1, w2, w3, u1, u2, u3;
 +      float v1[3], v2[3], relativeVelocity[3];
 +      float magrelVel;
 +      
 +      cloth1 = clmd->clothObject;
 +      cloth2 = coll_clmd->clothObject;
 +
 +      // search = clmd->coll_parms.collision_list;
 +      
 +      while(search)
 +      {
 +      collpair = search->link;
 +              
 +              // compute barycentric coordinates for both collision points
 +      collisions_compute_barycentric(collpair->pa,
 +      cloth1->verts[collpair->ap1].txold,
 +      cloth1->verts[collpair->ap2].txold,
 +      cloth1->verts[collpair->ap3].txold, 
 +      &w1, &w2, &w3);
 +      
 +      collisions_compute_barycentric(collpair->pb,
 +      cloth2->verts[collpair->bp1].txold,
 +      cloth2->verts[collpair->bp2].txold,
 +      cloth2->verts[collpair->bp3].txold,
 +      &u1, &u2, &u3);
 +      
 +              // Calculate relative "velocity".
 +      interpolateOnTriangle(v1, cloth1->verts[collpair->ap1].tv, cloth1->verts[collpair->ap2].tv, cloth1->verts[collpair->ap3].tv, w1, w2, w3);
 +              
 +      interpolateOnTriangle(v2, cloth2->verts[collpair->bp1].tv, cloth2->verts[collpair->bp2].tv, cloth2->verts[collpair->bp3].tv, u1, u2, u3);
 +              
 +      VECSUB(relativeVelocity, v1, v2);
 +                      
 +              // Calculate the normal component of the relative velocity (actually only the magnitude - the direction is stored in 'normal').
 +      magrelVel = INPR(relativeVelocity, collpair->normal);
 +              
 +              // printf("magrelVel: %f\n", magrelVel);
 +                              
 +              // Calculate masses of points.
 +              
 +              // If v_n_mag < 0 the edges are approaching each other.
 +      if(magrelVel < -ALMOST_ZERO) 
 +      {
 +                      // Calculate Impulse magnitude to stop all motion in normal direction.
 +                      // const double I_mag = v_n_mag / (1/m1 + 1/m2);
 +      float magnitude_i = magrelVel / 2.0f; // TODO implement masses
 +      float tangential[3], magtangent, magnormal, collvel[3];
 +      float vrel_t_pre[3];
 +      float vrel_t[3];
 +      double impulse;
 +      float epsilon = clmd->coll_parms.epsilon;
 +      float overlap = (epsilon + ALMOST_ZERO-collpair->distance);
 +                      
 +                      // calculateFrictionImpulse(tangential, relativeVelocity, collpair->normal, magrelVel, clmd->coll_parms.friction*0.01, magrelVel);
 +                      
 +                      // magtangent = INPR(tangential, tangential);
 +                      
 +                      // Apply friction impulse.
 +      if (magtangent < -ALMOST_ZERO) 
 +      {
 +                              
 +                              // printf("friction applied: %f\n", magtangent);
 +                              // TODO check original code 
 +}
 +                      
 +
 +      impulse = -2.0f * magrelVel / ( 1.0 + w1*w1 + w2*w2 + w3*w3);
 +                      
 +                      // printf("impulse: %f\n", impulse);
 +                      
 +                      // face A
 +      VECADDMUL(cloth1->verts[collpair->ap1].impulse, collpair->normal, w1 * impulse); 
 +      cloth1->verts[collpair->ap1].impulse_count++;
 +                      
 +      VECADDMUL(cloth1->verts[collpair->ap2].impulse, collpair->normal, w2 * impulse); 
 +      cloth1->verts[collpair->ap2].impulse_count++;
 +                      
 +      VECADDMUL(cloth1->verts[collpair->ap3].impulse, collpair->normal, w3 * impulse); 
 +      cloth1->verts[collpair->ap3].impulse_count++;
 +                      
 +                      // face B
 +      VECADDMUL(cloth2->verts[collpair->bp1].impulse, collpair->normal, u1 * impulse); 
 +      cloth2->verts[collpair->bp1].impulse_count++;
 +                      
 +      VECADDMUL(cloth2->verts[collpair->bp2].impulse, collpair->normal, u2 * impulse); 
 +      cloth2->verts[collpair->bp2].impulse_count++;
 +                      
 +      VECADDMUL(cloth2->verts[collpair->bp3].impulse, collpair->normal, u3 * impulse); 
 +      cloth2->verts[collpair->bp3].impulse_count++;
 +                      
 +                      
 +      result = 1;     
 +              
 +                      // printf("magnitude_i: %f\n", magnitude_i); // negative before collision in my case
 +                      
 +                      // Apply the impulse and increase impulse counters.
 +      
 +                      
 +}
 +              
 +      search = search->next;
 +}
 +      
 +      
 +      return result;
 +      */
 +      return 0;
 +}
 +
 +
 +int collisions_collision_response_moving_tris(ClothModifierData *clmd, ClothModifierData *coll_clmd)
 +{
 +      return 0;
 +}
 +
 +
 +int collisions_collision_response_moving_edges(ClothModifierData *clmd, ClothModifierData *coll_clmd)
 +{
 +      return 0;
 +}
 +
 +void cloth_collision_static(ClothModifierData *clmd, LinkNode *collision_list)
 +{
 +      /*
 +      CollPair *collpair = NULL;
 +      Cloth *cloth1=NULL, *cloth2=NULL;
 +      MFace *face1=NULL, *face2=NULL;
 +      ClothVertex *verts1=NULL, *verts2=NULL;
 +      double distance = 0;
 +      float epsilon = clmd->coll_parms->epsilon;
 +      unsigned int i = 0;
 +
 +      for(i = 0; i < 4; i++)
 +      {
 +      collpair = (CollPair *)MEM_callocN(sizeof(CollPair), "cloth coll pair");        
 +              
 +      cloth1 = clmd->clothObject;
 +      cloth2 = coll_clmd->clothObject;
 +              
 +      verts1 = cloth1->verts;
 +      verts2 = cloth2->verts;
 +      
 +      face1 = &(cloth1->mfaces[tree1->tri_index]);
 +      face2 = &(cloth2->mfaces[tree2->tri_index]);
 +              
 +              // check all possible pairs of triangles
 +      if(i == 0)
 +      {
 +      collpair->ap1 = face1->v1;
 +      collpair->ap2 = face1->v2;
 +      collpair->ap3 = face1->v3;
 +                      
 +      collpair->bp1 = face2->v1;
 +      collpair->bp2 = face2->v2;
 +      collpair->bp3 = face2->v3;
 +                      
 +}
 +              
 +      if(i == 1)
 +      {
 +      if(face1->v4)
 +      {
 +      collpair->ap1 = face1->v3;
 +      collpair->ap2 = face1->v4;
 +      collpair->ap3 = face1->v1;
 +                              
 +      collpair->bp1 = face2->v1;
 +      collpair->bp2 = face2->v2;
 +      collpair->bp3 = face2->v3;
 +}
 +      else
 +      i++;
 +}
 +              
 +      if(i == 2)
 +      {
 +      if(face2->v4)
 +      {
 +      collpair->ap1 = face1->v1;
 +      collpair->ap2 = face1->v2;
 +      collpair->ap3 = face1->v3;
 +                              
 +      collpair->bp1 = face2->v3;
 +      collpair->bp2 = face2->v4;
 +      collpair->bp3 = face2->v1;
 +}
 +      else
 +      i+=2;
 +}
 +              
 +      if(i == 3)
 +      {
 +      if((face1->v4)&&(face2->v4))
 +      {
 +      collpair->ap1 = face1->v3;
 +      collpair->ap2 = face1->v4;
 +      collpair->ap3 = face1->v1;
 +                              
 +      collpair->bp1 = face2->v3;
 +      collpair->bp2 = face2->v4;
 +      collpair->bp3 = face2->v1;
 +}
 +      else
 +      i++;
 +}
 +      
 +              
 +      if(i < 4)
 +      {
 +                      // calc distance + normal       
 +      distance = plNearestPoints(
 +      verts1[collpair->ap1].txold, verts1[collpair->ap2].txold, verts1[collpair->ap3].txold, verts2[collpair->bp1].txold, verts2[collpair->bp2].txold, verts2[collpair->bp3].txold, collpair->pa,collpair->pb,collpair->vector);
 +                      
 +      if (distance <= (epsilon + ALMOST_ZERO))
 +      {
 +                              // printf("dist: %f\n", (float)distance);
 +                              
 +                              // collpair->face1 = tree1->tri_index;
 +                              // collpair->face2 = tree2->tri_index;
 +                              
 +                              // VECCOPY(collpair->normal, collpair->vector);
 +                              // Normalize(collpair->normal);
 +                              
 +                              // collpair->distance = distance;
 +                              
 +}
 +      else
 +      {
 +      MEM_freeN(collpair);
 +}
 +}
 +      else
 +      {
 +      MEM_freeN(collpair);
 +}
 +}
 +      */
 +}
 +
 +int collisions_are_edges_adjacent(ClothModifierData *clmd, ClothModifierData *coll_clmd, EdgeCollPair *edgecollpair)
 +{
 +      Cloth *cloth1, *cloth2;
 +      ClothVertex *verts1, *verts2;
 +      float temp[3];
 +       /*
 +      cloth1 = clmd->clothObject;
 +      cloth2 = coll_clmd->clothObject;
 +      
 +      verts1 = cloth1->verts;
 +      verts2 = cloth2->verts;
 +      
 +      VECSUB(temp, verts1[edgecollpair->p11].xold, verts2[edgecollpair->p21].xold);
 +      if(ABS(INPR(temp, temp)) < ALMOST_ZERO)
 +      return 1;
 +      
 +      VECSUB(temp, verts1[edgecollpair->p11].xold, verts2[edgecollpair->p22].xold);
 +      if(ABS(INPR(temp, temp)) < ALMOST_ZERO)
 +      return 1;
 +      
 +      VECSUB(temp, verts1[edgecollpair->p12].xold, verts2[edgecollpair->p21].xold);
 +      if(ABS(INPR(temp, temp)) < ALMOST_ZERO)
 +      return 1;
 +      
 +      VECSUB(temp, verts1[edgecollpair->p12].xold, verts2[edgecollpair->p22].xold);
 +      if(ABS(INPR(temp, temp)) < ALMOST_ZERO)
 +      return 1;
 +       */
 +      return 0;
 +}
 +
 +
 +void collisions_collision_moving_edges(ClothModifierData *clmd, ClothModifierData *coll_clmd, CollisionTree *tree1, CollisionTree *tree2)
 +{
 +      /*
 +      EdgeCollPair edgecollpair;
 +      Cloth *cloth1=NULL, *cloth2=NULL;
 +      MFace *face1=NULL, *face2=NULL;
 +      ClothVertex *verts1=NULL, *verts2=NULL;
 +      double distance = 0;
 +      float epsilon = clmd->coll_parms->epsilon;
 +      unsigned int i = 0, j = 0, k = 0;
 +      int numsolutions = 0;
 +      float a[3], b[3], c[3], d[3], e[3], f[3], solution[3];
 +      
 +      cloth1 = clmd->clothObject;
 +      cloth2 = coll_clmd->clothObject;
 +      
 +      verts1 = cloth1->verts;
 +      verts2 = cloth2->verts;
 +
 +      face1 = &(cloth1->mfaces[tree1->tri_index]);
 +      face2 = &(cloth2->mfaces[tree2->tri_index]);
 +      
 +      for( i = 0; i < 5; i++)
 +      {
 +      if(i == 0) 
 +      {
 +      edgecollpair.p11 = face1->v1;
 +      edgecollpair.p12 = face1->v2;
 +}
 +      else if(i == 1) 
 +      {
 +      edgecollpair.p11 = face1->v2;
 +      edgecollpair.p12 = face1->v3;
 +}
 +      else if(i == 2) 
 +      {
 +      if(face1->v4) 
 +      {
 +      edgecollpair.p11 = face1->v3;
 +      edgecollpair.p12 = face1->v4;
 +}
 +      else 
 +      {
 +      edgecollpair.p11 = face1->v3;
 +      edgecollpair.p12 = face1->v1;
 +      i+=5; // get out of here after this edge pair is handled
 +}
 +}
 +      else if(i == 3) 
 +      {
 +      if(face1->v4) 
 +      {
 +      edgecollpair.p11 = face1->v4;
 +      edgecollpair.p12 = face1->v1;
 +}     
 +      else
 +      continue;
 +}
 +      else
 +      {
 +      edgecollpair.p11 = face1->v3;
 +      edgecollpair.p12 = face1->v1;
 +}
 +
 +              
 +      for( j = 0; j < 5; j++)
 +      {
 +      if(j == 0)
 +      {
 +      edgecollpair.p21 = face2->v1;
 +      edgecollpair.p22 = face2->v2;
 +}
 +      else if(j == 1)
 +      {
 +      edgecollpair.p21 = face2->v2;
 +      edgecollpair.p22 = face2->v3;
 +}
 +      else if(j == 2)
 +      {
 +      if(face2->v4) 
 +      {
 +      edgecollpair.p21 = face2->v3;
 +      edgecollpair.p22 = face2->v4;
 +}
 +      else 
 +      {
 +      edgecollpair.p21 = face2->v3;
 +      edgecollpair.p22 = face2->v1;
 +}
 +}
 +      else if(j == 3)
 +      {
 +      if(face2->v4) 
 +      {
 +      edgecollpair.p21 = face2->v4;
 +      edgecollpair.p22 = face2->v1;
 +}
 +      else
 +      continue;
 +}
 +      else
 +      {
 +      edgecollpair.p21 = face2->v3;
 +      edgecollpair.p22 = face2->v1;
 +}
 +                      
 +                      
 +      if(!collisions_are_edges_adjacent(clmd, coll_clmd, &edgecollpair))
 +      {
 +      VECSUB(a, verts1[edgecollpair.p12].xold, verts1[edgecollpair.p11].xold);
 +      VECSUB(b, verts1[edgecollpair.p12].v, verts1[edgecollpair.p11].v);
 +      VECSUB(c, verts1[edgecollpair.p21].xold, verts1[edgecollpair.p11].xold);
 +      VECSUB(d, verts1[edgecollpair.p21].v, verts1[edgecollpair.p11].v);
 +      VECSUB(e, verts2[edgecollpair.p22].xold, verts1[edgecollpair.p11].xold);
 +      VECSUB(f, verts2[edgecollpair.p22].v, verts1[edgecollpair.p11].v);
 +                              
 +      numsolutions = collisions_get_collision_time(a, b, c, d, e, f, solution);
 +                              
 +      for (k = 0; k < numsolutions; k++) 
 +      {                                                               
 +      if ((solution[k] >= 0.0) && (solution[k] <= 1.0)) 
 +      {
 +      float out_collisionTime = solution[k];
 +                                              
 +                                              // TODO: check for collisions 
 +                                              
 +                                              // TODO: put into (edge) collision list
 +                                              
 +      printf("Moving edge found!\n");
 +}
 +}
 +}
 +}
 +}     
 +      */      
 +}
 +
 +void collisions_collision_moving_tris(ClothModifierData *clmd, ClothModifierData *coll_clmd, CollisionTree *tree1, CollisionTree *tree2)
 +{
 +      /*
 +      CollPair collpair;
 +      Cloth *cloth1=NULL, *cloth2=NULL;
 +      MFace *face1=NULL, *face2=NULL;
 +      ClothVertex *verts1=NULL, *verts2=NULL;
 +      double distance = 0;
 +      float epsilon = clmd->coll_parms->epsilon;
 +      unsigned int i = 0, j = 0, k = 0;
 +      int numsolutions = 0;
 +      float a[3], b[3], c[3], d[3], e[3], f[3], solution[3];
 +
 +      for(i = 0; i < 2; i++)
 +      {               
 +      cloth1 = clmd->clothObject;
 +      cloth2 = coll_clmd->clothObject;
 +              
 +      verts1 = cloth1->verts;
 +      verts2 = cloth2->verts;
 +      
 +      face1 = &(cloth1->mfaces[tree1->tri_index]);
 +      face2 = &(cloth2->mfaces[tree2->tri_index]);
 +              
 +              // check all possible pairs of triangles
 +      if(i == 0)
 +      {
 +      collpair.ap1 = face1->v1;
 +      collpair.ap2 = face1->v2;
 +      collpair.ap3 = face1->v3;
 +                      
 +      collpair.pointsb[0] = face2->v1;
 +      collpair.pointsb[1] = face2->v2;
 +      collpair.pointsb[2] = face2->v3;
 +      collpair.pointsb[3] = face2->v4;
 +}
 +              
 +      if(i == 1)
 +      {
 +      if(face1->v4)
 +      {
 +      collpair.ap1 = face1->v3;
 +      collpair.ap2 = face1->v4;
 +      collpair.ap3 = face1->v1;
 +                              
 +      collpair.pointsb[0] = face2->v1;
 +      collpair.pointsb[1] = face2->v2;
 +      collpair.pointsb[2] = face2->v3;
 +      collpair.pointsb[3] = face2->v4;
 +}
 +      else
 +      i++;
 +}
 +              
 +              // calc SIPcode (?)
 +              
 +      if(i < 2)
 +      {
 +      VECSUB(a, verts1[collpair.ap2].xold, verts1[collpair.ap1].xold);
 +      VECSUB(b, verts1[collpair.ap2].v, verts1[collpair.ap1].v);
 +      VECSUB(c, verts1[collpair.ap3].xold, verts1[collpair.ap1].xold);
 +      VECSUB(d, verts1[collpair.ap3].v, verts1[collpair.ap1].v);
 +                              
 +      for(j = 0; j < 4; j++)
 +      {                                       
 +      if((j==3) && !(face2->v4))
 +      break;
 +                              
 +      VECSUB(e, verts2[collpair.pointsb[j]].xold, verts1[collpair.ap1].xold);
 +      VECSUB(f, verts2[collpair.pointsb[j]].v, verts1[collpair.ap1].v);
 +                              
 +      numsolutions = collisions_get_collision_time(a, b, c, d, e, f, solution);
 +                              
 +      for (k = 0; k < numsolutions; k++) 
 +      {                                                               
 +      if ((solution[k] >= 0.0) && (solution[k] <= 1.0)) 
 +      {
 +      float out_collisionTime = solution[k];
 +                                              
 +                                              // TODO: check for collisions 
 +                                              
 +                                              // TODO: put into (point-face) collision list
 +                                              
 +      printf("Moving found!\n");
 +                                              
 +}
 +}
 +                              
 +                              // TODO: check borders for collisions
 +}
 +                      
 +}
 +}
 +      */
 +}
 +
 +void collisions_collision_moving(ClothModifierData *clmd, ClothModifierData *coll_clmd, CollisionTree *tree1, CollisionTree *tree2)
 +{
 +      /*
 +      // TODO: check for adjacent
 +      collisions_collision_moving_edges(clmd, coll_clmd, tree1, tree2);
 +      
 +      collisions_collision_moving_tris(clmd, coll_clmd, tree1, tree2);
 +      collisions_collision_moving_tris(coll_clmd, clmd, tree2, tree1);
 +      */
 +}
 +
 +// cloth - object collisions
 +int cloth_bvh_objcollision(ClothModifierData * clmd, float step, float prevstep, float dt)
 +{
 +      
 +      Base *base = NULL;
 +      CollisionModifierData *collmd = NULL;
 +      Cloth *cloth = NULL;
 +      Object *ob2 = NULL;
 +      BVH *bvh1 = NULL, *bvh2 = NULL, *self_bvh;
 +      LinkNode *collision_list = NULL; 
 +      unsigned int i = 0, j = 0;
 +      int collisions = 0, count = 0;
 +      float (*current_x)[3];
 +
 +      if (!(((Cloth *)clmd->clothObject)->tree))
 +      {
 +              printf("No BVH found\n");
 +              return 0;
 +      }
 +      
 +      cloth = clmd->clothObject;
 +      bvh1 = cloth->tree;
 +      self_bvh = cloth->selftree;
 +      
 +      ////////////////////////////////////////////////////////////
 +      // static collisions
 +      ////////////////////////////////////////////////////////////
 +      
 +      // update cloth bvh
 +      bvh_update_from_float3(bvh1, cloth->current_xold, cloth->numverts, cloth->current_x, 0); // 0 means STATIC, 1 means MOVING (see later in this function)
 +/*
 +      // check all collision objects
 +      for (base = G.scene->base.first; base; base = base->next)
 +      {
 +      ob2 = base->object;
 +      collmd = (CollisionModifierData *) modifiers_findByType (ob2, eModifierType_Collision);
 +              
 +      if (!collmd)
 +      continue;
 +              
 +              // check if there is a bounding volume hierarchy
 +      if (collmd->tree) 
 +      {                       
 +      bvh2 = collmd->tree;
 +                      
 +                      // update position + bvh of collision object
 +      collision_move_object(collmd, step, prevstep);
 +      bvh_update_from_mvert(collmd->tree, collmd->current_x, collmd->numverts, NULL, 0);
 +                      
 +                      // fill collision list 
 +      collisions += bvh_traverse(bvh1->root, bvh2->root, &collision_list);
 +                      
 +                      // call static collision response
 +                      
 +                      // free collision list
 +      if(collision_list)
 +      {
 +      LinkNode *search = collision_list;
 +                              
 +      while(search)
 +      {
 +      CollisionPair *coll_pair = search->link;
 +                                      
 +      MEM_freeN(coll_pair);
 +      search = search->next;
 +}
 +      BLI_linklist_free(collision_list,NULL);
 +      
 +      collision_list = NULL;
 +}
 +}
 +}
 +      
 +      //////////////////////////////////////////////
 +      // update velocities + positions
 +      //////////////////////////////////////////////
 +      for(i = 0; i < cloth->numverts; i++)
 +      {
 +      VECADD(cloth->current_x[i], cloth->current_xold[i], cloth->current_v[i]);
 +}
 +      //////////////////////////////////////////////
 +*/    
 +      /*
 +      // fill collision list 
 +      collisions += bvh_traverse(self_bvh->root, self_bvh->root, &collision_list);
 +      
 +      // call static collision response
 +      
 +      // free collision list
 +      if(collision_list)
 +      {
 +      LinkNode *search = collision_list;
 +              
 +      while(search)
 +      {
 +      float distance = 0;
 +      float mindistance = cloth->selftree->epsilon;
 +      CollisionPair *collpair = (CollisionPair *)search->link;
 +                      
 +                      // get distance of faces
 +      distance = plNearestPoints(
 +      cloth->current_x[collpair->point_indexA[0]], cloth->current_x[collpair->point_indexA[1]], cloth->current_x[collpair->point_indexA[2]], cloth->current_x[collpair->point_indexB[0]], cloth->current_x[collpair->point_indexB[1]], cloth->current_x[collpair->point_indexB[2]], collpair->pa,collpair->pb,collpair->vector);
 +                                      
 +      if(distance < mindistance)
 +      {
 +      ///////////////////////////////////////////
 +                              // TODO: take velocity of the collision points into account!
 +      ///////////////////////////////////////////
 +                              
 +      float correction = mindistance - distance;
 +      float temp[3];
 +                              
 +      VECCOPY(temp, collpair->vector);
 +      Normalize(temp);
 +      VecMulf(temp, -correction*0.5);
 +                              
 +      if(!((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [collpair->point_indexA[0]].goal >= SOFTGOALSNAP)))
 +      VECSUB(cloth->current_x[collpair->point_indexA[0]], cloth->current_x[collpair->point_indexA[0]], temp); 
 +                              
 +      if(!((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [collpair->point_indexA[1]].goal >= SOFTGOALSNAP)))
 +      VECSUB(cloth->current_x[collpair->point_indexA[1]], cloth->current_x[collpair->point_indexA[1]], temp);
 +                              
 +      if(!((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [collpair->point_indexA[2]].goal >= SOFTGOALSNAP)))
 +      VECSUB(cloth->current_x[collpair->point_indexA[2]], cloth->current_x[collpair->point_indexA[2]], temp);
 +                              
 +                              
 +      if(!((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [collpair->point_indexB[0]].goal >= SOFTGOALSNAP)))
 +      VECSUB(cloth->current_x[collpair->point_indexB[0]], cloth->current_x[collpair->point_indexB[0]], temp); 
 +                              
 +      if(!((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [collpair->point_indexB[1]].goal >= SOFTGOALSNAP)))
 +      VECSUB(cloth->current_x[collpair->point_indexB[1]], cloth->current_x[collpair->point_indexB[1]], temp);
 +                              
 +      if(!((clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL) && (cloth->verts [collpair->point_indexB[2]].goal >= SOFTGOALSNAP)))
 +      VECSUB(cloth->current_x[collpair->point_indexB[2]], cloth->current_x[collpair->point_indexB[2]], temp);
 +                                      
 +      collisions = 1;
 +                              
 +}
 +                      
 +}
 +              
 +      search = collision_list;
 +      while(search)
 +      {
 +      CollisionPair *coll_pair = search->link;
 +                      
 +      MEM_freeN(coll_pair);
 +      search = search->next;
 +}
 +      BLI_linklist_free(collision_list,NULL);
 +
 +      collision_list = NULL;
 +}
 +      */
 +      // Test on *simple* selfcollisions
 +      collisions = 1;
 +      count = 0;
 +      current_x = cloth->current_x; // needed for openMP
-       VecMulf(temp, -correction*0.5);
-       VECADD(current_x[j], current_x[j], temp);
-                                               
-       VECSUB(current_x[i], current_x[i], temp);       
- }
++      
++      // #pragma omp parallel for private(i,j, collisions) shared(current_x)
++      // for ( count = 0; count < 6; count++ )
 +      {
-       collisions = 1;
- }
- }
- }
- }
++              collisions = 0;
++      
++              for ( i = 0; i < cloth->numverts; i++ )
++              {
++                      float mindistance1 = cloth->verts[i].collball;
++                      
++                      for ( j = i + 1; j < cloth->numverts; j++ )
++                      {
++                              float temp[3];
++                              float length = 0;
++                              
++                              float mindistance2 = cloth->verts[j].collball;
++      
++                              if ( clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL )
++                              {
++                                      if ( ( cloth->verts [i].goal >= SOFTGOALSNAP )
++                                              && ( cloth->verts [j].goal >= SOFTGOALSNAP ) )
++                                      {
++                                              continue;
++                                      }
++                              }
++      
++                              // check for adjacent points
++                              if ( BLI_edgehash_haskey ( cloth->edgehash, i, j ) )
++                              {
++                                      continue;
++                              }
++      
++                              VECSUB ( temp, current_x[i], current_x[j] );
++      
++                              length = Normalize ( temp );
++      
++                              if ( length < ((mindistance1 + mindistance2)) )
++                              {
++                                      float correction = ((mindistance1 + mindistance2)) - length;
 +                                      
-       for(i = 0; i < cloth->numverts; i++)
++                                      printf("correction: %f\n", correction);
++      
++                                      if ( ( clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL ) && ( cloth->verts [i].goal >= SOFTGOALSNAP ) )
++                                      {
++                                              VecMulf ( temp, -correction );
++                                              VECADD ( current_x[j], current_x[j], temp );
++                                      }
++                                      else if ( ( clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL ) && ( cloth->verts [j].goal >= SOFTGOALSNAP ) )
++                                      {
++                                              VecMulf ( temp, correction );
++                                              VECADD ( current_x[i], current_x[i], temp );
++                                      }
++                                      else
++                                      {
++                                              VecMulf ( temp, -correction*0.5 );
++                                              VECADD ( current_x[j], current_x[j], temp );
++      
++                                              VECSUB ( current_x[i], current_x[i], temp );
++                                      }
++      
++                                      collisions = 1;
++                              }
++                      }
++              }
++      }
++      
 +      
 +      //////////////////////////////////////////////
 +      // SELFCOLLISIONS: update velocities
 +      //////////////////////////////////////////////
-       VECSUB(cloth->current_v[i], cloth->current_x[i], cloth->current_xold[i]);
- }
++      for ( i = 0; i < cloth->numverts; i++ )
 +      {
- */    
++              VECSUB ( cloth->current_v[i], cloth->current_x[i], cloth->current_xold[i] );
++      }
 +      //////////////////////////////////////////////
++
 +      return 1;
 +}
index bd9d1cb75ca4be6058dbe5578bf266e3b883deae,bd9d1cb75ca4be6058dbe5578bf266e3b883deae..4a3d8f70b1df1c343fe8cf3f141898213e5d4456
@@@ -778,7 -778,7 +778,7 @@@ static void calculate_collision_balls(O
                        bs = sb->bspring + bp->springs[b-1];
                        if (bs->order == 1){
                        akku += bs->len;
--                      akku_count++,
++                      akku_count++;
                        min = MIN2(bs->len,min);
                        max = MAX2(bs->len,max);
                        }
index 6ae2cd73df3fca177f2df260e84873f8c2f6ca54,f356c97c987fdba83d0ac8d34b73f9e8e70b4ce3..959fd77095961c9f1ff10706b986ca7000098d3c
@@@ -1654,10 -1660,12 +1670,12 @@@ static void draw_modifier(uiBlock *bloc
                uiBlockBeginAlign(block);
                uiDefBut(block, TEX, B_MODIFIER_REDRAW, "", x+10, y-1, buttonWidth-60, 19, md->name, 0.0, sizeof(md->name)-1, 0.0, 0.0, "Modifier name"); 
  
 -                      /* Softbody not allowed in this situation, enforce! */
 -              if (md->type!=eModifierType_Softbody || !(ob->pd && ob->pd->deflect)) {
 +              /* Softbody not allowed in this situation, enforce! */
 +              if ((md->type!=eModifierType_Softbody && md->type!=eModifierType_Collision) || !(ob->pd && ob->pd->deflect)) {
                        uiDefIconButBitI(block, TOG, eModifierMode_Render, B_MODIFIER_RECALC, ICON_SCENE, x+10+buttonWidth-60, y-1, 19, 19,&md->mode, 0, 0, 1, 0, "Enable modifier during rendering");
-                       uiDefIconButBitI(block, TOG, eModifierMode_Realtime, B_MODIFIER_RECALC, VICON_VIEW3D, x+10+buttonWidth-40, y-1, 19, 19,&md->mode, 0, 0, 1, 0, "Enable modifier during interactive display");
+                       but= uiDefIconButBitI(block, TOG, eModifierMode_Realtime, B_MODIFIER_RECALC, VICON_VIEW3D, x+10+buttonWidth-40, y-1, 19, 19,&md->mode, 0, 0, 1, 0, "Enable modifier during interactive display");
+                       if (md->type==eModifierType_ParticleSystem)
+                               uiButSetFunc(but, modifiers_psysEnable, ob, md);
                        if (mti->flags&eModifierTypeFlag_SupportsEditmode) {
                                uiDefIconButBitI(block, TOG, eModifierMode_Editmode, B_MODIFIER_RECALC, VICON_EDIT, x+10+buttonWidth-20, y-1, 19, 19,&md->mode, 0, 0, 1, 0, "Enable modifier during Editmode (only if enabled for display)");
                        }
index a5be0e0357f080ce5459730977eda8b2698d67cd,731fe49e4c9d945e68b32623ec4975e0e6a79794..2ce543eec443b17a534b00ff8cca0ff09241197a
@@@ -4620,233 -4618,6 +4620,233 @@@ errMessage
  #endif // DISABLE_ELBEEM
  }
  
-                       uiDefButF(block, NUM, B_CLOTH_RENEW, "Selfcoll balls:",    10,120,150,20, &clmd->coll_parms->selfepsilon, 0.001f, 1.0, 0.01f, 0, "Minimum distance between two selfcollision points");
 +/* Panel for cloth */
 +static void object_cloth__enabletoggle(void *ob_v, void *arg2)
 +{
 +      Object *ob = ob_v;
 +      ModifierData *md = modifiers_findByType(ob, eModifierType_Cloth);
 +
 +      if (!md) {
 +              md = modifier_new(eModifierType_Cloth);
 +              BLI_addhead(&ob->modifiers, md);
 +      }
 +      else {
 +              BLI_remlink(&ob->modifiers, md);
 +              modifier_free(md);
 +      }
 +
 +      allqueue(REDRAWBUTSEDIT, 0);
 +}
 +
 +static void object_panel_cloth(Object *ob)
 +{
 +      uiBlock *block;
 +      static int val, val2;
 +      uiBut *but;
 +      ClothModifierData *clmd = (ClothModifierData *)modifiers_findByType(ob, eModifierType_Cloth);
 +      block= uiNewBlock(&curarea->uiblocks, "object_panel_cloth", UI_EMBOSS, UI_HELV, curarea->win);
 +      if(uiNewPanel(curarea, block, "Cloth", "Physics", 640, 0, 318, 204)==0) return;
 +
 +      if(ob->id.lib) uiSetButLock(1, "Can't edit library data");
 +      val = ((clmd)?(1):(0));
 +
 +      but = uiDefButI(block, TOG, REDRAWBUTSOBJECT, "Cloth Object",   10,200,130,20, &val, 0, 0, 0, 0, "Sets object to become cloth");
 +      uiButSetFunc(but, object_cloth__enabletoggle, ob, NULL);
 +      uiDefBut(block, LABEL, 0, "",10,10,300,0, NULL, 0.0, 0, 0, 0, ""); /* tell UI we go to 10,10*/
 +      
 +      if(clmd)
 +      {
 +              // but = uiDefButBitI(block, TOG, CLOTH_SIMSETTINGS_FLAG_COLLOBJ, B_EFFECT_DEP, "Collision Object",     170,200,130,20, &clmd->sim_parms->flags, 0, 0, 0, 0, "Sets object to become a cloth collision object");
 +
 +              if (!(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_COLLOBJ))
 +              {
 +                      Cloth *cloth = clmd->clothObject;
 +                      int defCount;
 +                      char *clvg1, *clvg2;
 +                      char clmvg [] = "Mass Vertex Group%t|None%x0|";
 +      
 +                      val2=0;
 +      
 +      //              uiDefButBitI(block, TOG, CSIMSETT_FLAG_ADVANCED, REDRAWBUTSOBJECT, "Advanced",  180,200,130,20, &clmd->sim_parms->flags, 0, 0, 0, 0, "Enable advanced mode");
 +                      
 +                      /* GENERAL STUFF */
 +                      uiClearButLock();
 +                      uiBlockBeginAlign(block);
 +                      uiDefButF(block, NUM, B_CLOTH_RENEW, "StructStiff:",       10,170,150,20, &clmd->sim_parms->structural, 1.0, 10000.0, 100, 0, "Overall stiffness of structure");
 +                      uiDefButF(block, NUM, B_CLOTH_RENEW, "BendStiff:",         160,170,150,20, &clmd->sim_parms->bending, 0.0, 10000.0, 1000, 0, "Wrinkle possibility");
 +                      uiDefButI(block, NUM, B_CLOTH_RENEW, "Quality:",  10,150,150,20, &clmd->sim_parms->stepsPerFrame, 3.0, 10.0, 5, 0, "Quality of the simulation (higher=>better=>slower)");
 +                      uiBlockEndAlign(block);
 +                      uiBlockBeginAlign(block);
 +                      uiDefButF(block, NUM, B_CLOTH_RENEW, "Spring Damp:",       160,150,150,20, &clmd->sim_parms->Cdis, 0.0, 10.0, 10, 0, "Spring damping");
 +                      uiDefButF(block, NUM, B_DIFF, "Air Damp:",         10,130,150,20, &clmd->sim_parms->Cvi, 0.0, 10.0, 10, 0, "Apply gravitation to point movement");
 +                      uiBlockEndAlign(block);                 
 +                      
 +                      uiClearButLock();
 +                      
 +                      uiBlockBeginAlign(block);
 +                      uiDefBut(block, LABEL, 0, "Gravity:",  10,100,60,20, NULL, 0.0, 0, 0, 0, "");
 +                      // uiClearButLock();
 +                      
 +                      uiDefButF(block, NUM, B_CLOTH_RENEW, "X:",         70,100,80,20, &clmd->sim_parms->gravity[0], -100.0, 100.0, 10, 0, "Apply gravitation to point movement");
 +                      uiDefButF(block, NUM, B_CLOTH_RENEW, "Y:",         150,100,80,20, &clmd->sim_parms->gravity[1], -100.0, 100.0, 10, 0, "Apply gravitation to point movement");
 +                      uiDefButF(block, NUM, B_CLOTH_RENEW, "Z:",         230,100,80,20, &clmd->sim_parms->gravity[2], -100.0, 100.0, 10, 0, "Apply gravitation to point movement");
 +                      uiBlockEndAlign(block);
 +                      
 +                      /* GOAL STUFF */
 +                      uiBlockBeginAlign(block);
 +                      uiDefButBitI(block, TOG, CLOTH_SIMSETTINGS_FLAG_GOAL, REDRAWVIEW3D, "Use Goal", 10,70,130,20, &clmd->sim_parms->flags, 0, 0, 0, 0, "Define forces for vertices to stick to animated position");
 +                      if (clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_GOAL)
 +                      {
 +                              if(ob->type==OB_MESH) 
 +                              {
 +                                      
 +                                      defCount = sizeof (clmvg);
 +                                      clvg1 = get_vertexgroup_menustr (ob);
 +                                      clvg2 = MEM_callocN (strlen (clvg1) + 1 + defCount, "clothVgMS");
 +                                      if (! clvg2) {
 +                                              printf ("draw_modifier: error allocating memory for cloth vertex group menu string.\n");
 +                                              return;
 +                                      }
 +                                      defCount = BLI_countlist (&ob->defbase);
 +                                      if (defCount == 0) 
 +                                      {
 +                                              clmd->sim_parms->vgroup_mass = 0;
 +                                      }
 +                                      sprintf (clvg2, "%s%s", clmvg, clvg1);
 +                                      
 +                                      uiDefButS(block, MENU, B_CLOTH_RENEW, clvg2,    140,70,20,20, &clmd->sim_parms->vgroup_mass, 0, defCount, 0, 0, "Browses available vertex groups");     
 +                                      MEM_freeN (clvg1);
 +                                      MEM_freeN (clvg2);
 +                                      
 +                                      if(clmd->sim_parms->vgroup_mass) 
 +                                      {
 +                                              bDeformGroup *defGroup = BLI_findlink(&ob->defbase, clmd->sim_parms->vgroup_mass-1);
 +                                              if(defGroup)
 +                                                      uiDefBut(block, BUT, B_DIFF, defGroup->name,    160,70,130,20, NULL, 0.0, 0.0, 0, 0, "Name of current vertex group");
 +                                              else
 +                                                      uiDefBut(block, BUT, B_DIFF, "(no group)",      160,70,130,20, NULL, 0.0, 0.0, 0, 0, "Vertex Group doesn't exist anymore");
 +                                              
 +                                              uiDefIconBut(block, BUT, B_CLOTH_DEL_VG, ICON_X, 290,70,20,20, 0, 0, 0, 0, 0, "Disable use of vertex group");
 +                                              
 +                                      }
 +                                      else
 +                                              uiDefButF(block, NUM, B_CLOTH_RENEW, "Goal:",   160,70,150,20, &clmd->sim_parms->defgoal, 0.0, 1.0, 10, 0, "Default Goal (vertex target position) value, when no Vertex Group used");
 +                              
 +                              }
 +                              else 
 +                              {
 +                                      uiDefButS(block, TOG, B_CLOTH_RENEW, "W",                       140,70,20,20, &clmd->sim_parms->vgroup_mass, 0, 1, 0, 0, "Use control point weight values");
 +                                      uiDefButF(block, NUM, B_CLOTH_RENEW, "Goal:",   160,70,150,20, &clmd->sim_parms->defgoal, 0.0, 1.0, 10, 0, "Default Goal (vertex target position) value, when no Vertex Group used");
 +                              }
 +                              
 +                              uiDefButF(block, NUM, B_CLOTH_RENEW, "G Stiff:",        10,50,150,20, &clmd->sim_parms->goalspring, 0.0, 500.0, 10, 0, "Goal (vertex target position) spring stiffness");
 +                              uiDefButF(block, NUM, B_CLOTH_RENEW, "G Damp:", 160,50,150,20, &clmd->sim_parms->goalfrict, 0.0, 50.0, 10, 0, "Goal (vertex target position) friction");
 +                              uiDefButF(block, NUM, B_CLOTH_RENEW, "G Min:",          10,30,150,20, &clmd->sim_parms->mingoal, 0.0, 1.0, 10, 0, "Goal minimum, vertex group weights are scaled to match this range");
 +                              uiDefButF(block, NUM, B_CLOTH_RENEW, "G Max:",          160,30,150,20, &clmd->sim_parms->maxgoal, 0.0, 1.0, 10, 0, "Goal maximum, vertex group weights are scaled to match this range");
 +                      }
 +                      uiBlockEndAlign(block); 
 +                      
 +                      /*
 +                      // no tearing supported anymore since modifier stack restrictions 
 +                      uiBlockBeginAlign(block);
 +                      uiDefButBitI(block, TOG, CSIMSETT_FLAG_TEARING_ENABLED, B_EFFECT_DEP, "Tearing",        10,0,150,20, &clmd->sim_parms->flags, 0, 0, 0, 0, "Sets object to become a cloth collision object");
 +                      
 +                      if (clmd->sim_parms->flags & CSIMSETT_FLAG_TEARING_ENABLED)
 +                      {
 +                              uiDefButI(block, NUM, B_DIFF, "Max extent:",       160,0,150,20, &clmd->sim_parms->maxspringlen, 1.0, 1000.0, 10, 0, "Maximum extension before spring gets cut");
 +                      }
 +                      
 +                      uiBlockEndAlign(block); 
 +                      */
 +              }
 +      }
 +}
 +
 +
 +static void object_panel_cloth_II(Object *ob)
 +{
 +      uiBlock *block;
 +      static int val;
 +      uiBut *but;
 +      ClothModifierData *clmd = NULL;
 +      
 +      clmd = (ClothModifierData *)modifiers_findByType(ob, eModifierType_Cloth);
 +      if(clmd)
 +      {
 +              if (!(clmd->sim_parms->flags & CLOTH_SIMSETTINGS_FLAG_COLLOBJ))
 +              {
 +                      Cloth *cloth = clmd->clothObject;
 +                      char str[128];
 +                      
 +                      block= uiNewBlock(&curarea->uiblocks, "object_panel_cloth_II", UI_EMBOSS, UI_HELV, curarea->win);
 +                      uiNewPanelTabbed("Cloth", "Physics");
 +                      if(uiNewPanel(curarea, block, "Cloth Cache", "Physics", 651, 0, 318, 204)==0) return;
 +              
 +                      uiSetButLock(object_data_is_libdata(ob), ERROR_LIBDATA_MESSAGE);
 +                      
 +                      uiDefButI(block, NUM, B_DIFF, "First Frame:",           10,160,150,20, &clmd->sim_parms->firstframe, 0, MAXFRAME, 1, 0, "Frame on which the simulation starts");
 +                      uiDefButI(block, NUM, B_DIFF, "Last Frame:",            160,160,150,20, &clmd->sim_parms->lastframe, 0, MAXFRAME, 10, 0, "Frame on which the simulation stops");
 +                      /*
 +                      if(clmd->sim_parms->cache)
 +                      {
 +                              int length = BLI_linklist_length(clmd->sim_parms->cache);
 +                              
 +                              // correct spelling if only 1 frame cacheed --> only gimmick  
 +                              if(length-clmd->sim_parms->preroll>1)
 +                                      sprintf (str, "Frame 1 - %d cached. [%d in preroll, %d in total]", length-clmd->sim_parms->preroll, clmd->sim_parms->preroll, length);
 +                              else
 +                                      sprintf (str, "Frame %d cached. [%d in preroll, %d in total]", length-clmd->sim_parms->preroll, clmd->sim_parms->preroll, length);
 +                              
 +                              uiDefBut(block, LABEL, 0, str,  10,140,290,20, NULL, 0.0, 0, 0, 0, "");
 +                              uiDefBut(block, LABEL, 0, "Clear cache:",  10,120,290,20, NULL, 0.0, 0, 0, 0, "");
 +                              uiBlockBeginAlign (block);
 +                              uiDefBut(block, BUT, B_CLOTH_CLEARCACHEALL, "All", 10, 100,145,20, NULL, 0.0, 0.0, 0, 0, "Free cloth cache without preroll");
 +                              uiDefBut(block, BUT, B_CLOTH_CLEARCACHEFRAME, "From next frame", 155, 100,145,20, NULL, 0.0, 0.0, 0, 0, "Free cloth cache");    
 +                              if(length>1) // B_CLOTH_CHANGEPREROLL
 +                                      uiDefButI(block, NUM, B_CLOTH_CHANGEPREROLL, "Preroll:", 10,80,145,20, &clmd->sim_parms->preroll, 0, length-1, 1, 0, "Simulation starts on this frame");        
 +                              else
 +                                      uiDefBut(block, LABEL, 0, " ",  10,80,145,20, NULL, 0.0, 0, 0, 0, "");
 +                      }
 +                      else 
 +                      {
 +                              uiDefBut(block, LABEL, 0, "No frames cached.",  10,120,290,20, NULL, 0.0, 0, 0, 0, "");
 +                      }*/
 +                      uiDefButBitI(block, TOG, CLOTH_SIMSETTINGS_FLAG_CCACHE_PROTECT, REDRAWVIEW3D, "Protect Cache",  10,50,145,20, &clmd->sim_parms->flags, 0, 0, 0, 0, "Protect cache from automatic freeing when scene changed");
 +                      uiBlockEndAlign(block);
 +              }
 +      }
 +      // uiBlockEndAlign(block);
 +}
 +
 +static void object_panel_cloth_III(Object *ob)
 +{
 +      uiBlock *block;
 +      ClothModifierData *clmd = NULL;
 +      
 +      clmd = (ClothModifierData *)modifiers_findByType(ob, eModifierType_Cloth);
 +      if(clmd)
 +      {       
 +              block= uiNewBlock(&curarea->uiblocks, "object_panel_cloth_III", UI_EMBOSS, UI_HELV, curarea->win);
 +              uiNewPanelTabbed("Cloth", "Physics");
 +              if(uiNewPanel(curarea, block, "Cloth Collisions", "Physics", 651, 0, 318, 204)==0) return;
 +      
 +              uiSetButLock(object_data_is_libdata(ob), ERROR_LIBDATA_MESSAGE);
 +              
 +              uiBlockBeginAlign(block);
 +              uiDefButBitI(block, TOG, CLOTH_COLLISIONSETTINGS_FLAG_ENABLED, REDRAWVIEW3D, "Enable collisions",       10,160,130,20, &clmd->coll_parms->flags, 0, 0, 0, 0, "Enable collisions with this object");
 +              if (clmd->coll_parms->flags & CLOTH_COLLISIONSETTINGS_FLAG_ENABLED)
 +              {
 +                      // uiDefBut(block, LABEL, 0, "",10,10,300,20, NULL, 0.0, 0, 0, 0, ""); /* tell UI we go to 10,10*/
 +                      uiDefButF(block, NUM, B_CLOTH_RENEW, "Min Distance:",      10,140,150,20, &clmd->coll_parms->epsilon, 0.001f, 1.0, 0.01f, 0, "Minimum distance between collision objects before collision response takes in");
 +                      uiDefBut(block, LABEL, 0, "",160,140,150,20, NULL, 0.0, 0, 0, 0, "");
++                      uiDefButF(block, NUM, B_CLOTH_RENEW, "Selfcoll balls:",    10,120,150,20, &clmd->coll_parms->selfepsilon, 0.001f, 1.0, 0.49f, 0, "Minimum distance between two selfcollision points");
 +              }
 +              else
 +                      uiDefBut(block, LABEL, 0, "",140,10,170,20, NULL, 0.0, 0, 0, 0, "");
 +              uiBlockEndAlign(block);
 +      }
 +}
 +
  void object_panels()
  {
        Object *ob;