diff options
author | Jens Ole Wund <bjornmose@gmx.net> | 2005-01-13 01:28:13 +0300 |
---|---|---|
committer | Jens Ole Wund <bjornmose@gmx.net> | 2005-01-13 01:28:13 +0300 |
commit | 912ef80bdc464de00fd2452c51bcd5cd3da050d6 (patch) | |
tree | 44f8e94182b52bf15ada6291d8d99839c4ce9af1 /source/blender/blenkernel | |
parent | ea740dea0251038978dc00f03b6b7c4c28e14314 (diff) |
big softbody commit
some vertex group , weight painting stuff too
/me crosses fingers it does not break anything
Diffstat (limited to 'source/blender/blenkernel')
-rw-r--r-- | source/blender/blenkernel/BKE_softbody.h | 6 | ||||
-rw-r--r-- | source/blender/blenkernel/intern/deform.c | 14 | ||||
-rw-r--r-- | source/blender/blenkernel/intern/displist.c | 2 | ||||
-rw-r--r-- | source/blender/blenkernel/intern/softbody.c | 738 |
4 files changed, 703 insertions, 57 deletions
diff --git a/source/blender/blenkernel/BKE_softbody.h b/source/blender/blenkernel/BKE_softbody.h index 77a68317818..24abff44a2b 100644 --- a/source/blender/blenkernel/BKE_softbody.h +++ b/source/blender/blenkernel/BKE_softbody.h @@ -35,8 +35,10 @@ #define BKE_SOFTBODY_H typedef struct BodyPoint { - float orig[3], pos[3], vec[3], force[3]; + float origS[3], origE[3], origT[3], pos[3], vec[3], force[3]; float weight, goal; + float prevpos[3], prevvec[3], prevdx[3], prevdv[3]; /* used for Heun integration */ + int nofsprings; int *springs; } BodyPoint; typedef struct BodySpring { @@ -57,7 +59,7 @@ typedef struct SoftBody { extern void free_softbody(SoftBody *sb); /* makes totally fresh start situation */ -extern void object_to_softbody(Object *ob); +extern void object_to_softbody(Object *ob,float ctime); /* copy original (but new) situation in softbody, as result of matrices or deform */ void object_update_softbody(Object *ob); diff --git a/source/blender/blenkernel/intern/deform.c b/source/blender/blenkernel/intern/deform.c index 16145ca2d63..19c2f15e7d3 100644 --- a/source/blender/blenkernel/intern/deform.c +++ b/source/blender/blenkernel/intern/deform.c @@ -227,7 +227,7 @@ int mesh_modifier(Object *ob, char mode) else if(ob->effect.first); // weak... particles too else if(ob->parent && ob->parent->type==OB_LATTICE); else if(ob->parent && ob->partype==PARSKEL); - else if(ob->softflag); + else if(ob->softflag & 0x01); else return 0; if(me->totvert==0) return 0; @@ -249,15 +249,21 @@ int mesh_modifier(Object *ob, char mode) } if(ob->effect.first) done |= object_wave(ob); - - if(ob->softflag) { + + if((ob->softflag & 0x01) && !(ob->softflag & 0x08)) { float ctime= bsystem_time(ob, NULL, (float)G.scene->r.cfra, 0.0); done= 1; object_softbody_step(ob, ctime); } - + /* deform: input mesh, output ob dl_verts. is used by subsurf (output should be in mesh ton!) */ done |= object_deform(ob); + + if((ob->softflag & 0x01) && (ob->softflag & 0x08)) { + float ctime= bsystem_time(ob, NULL, (float)G.scene->r.cfra, 0.0); + done= 1; + object_softbody_step(ob, ctime); + } /* put deformed vertices in dl->verts, optional subsurf will replace that */ if(done) { diff --git a/source/blender/blenkernel/intern/displist.c b/source/blender/blenkernel/intern/displist.c index 228588e2c4a..9c24e6361a0 100644 --- a/source/blender/blenkernel/intern/displist.c +++ b/source/blender/blenkernel/intern/displist.c @@ -2405,7 +2405,7 @@ void test_all_displists(void) } } - if(ob->softflag) freedisplist_object(ob); + if(ob->softflag & 0x01) freedisplist_object(ob); /* warn, ob pointer changed in case of OB_MALL */ if ELEM(ob->type, OB_CURVE, OB_SURF) { diff --git a/source/blender/blenkernel/intern/softbody.c b/source/blender/blenkernel/intern/softbody.c index 0c785041bb9..b56dec5417b 100644 --- a/source/blender/blenkernel/intern/softbody.c +++ b/source/blender/blenkernel/intern/softbody.c @@ -31,6 +31,30 @@ * ***** END GPL/BL DUAL LICENSE BLOCK ***** */ +/* +****** +variables on the UI for now +typedef struct Object { +..... + float formfactor, softtime; softtime = #euler integrations steps per frame +..... + float sb_goalspring; softbody goal springs + float sb_goalfrict; softbody goal springs friction + float sb_inspring; softbody inner springs + float sb_infrict; softbody inner springs friction + float sb_nodemass; softbody mass of *vertex* + float sb_grav; softbody amount of gravitaion to apply + float sb_mingoal; quick limits for goal + float sb_maxgoal; + float sb_mediafrict; friction to env + float sb_pad1; free + + + +***** +*/ + + #include <math.h> #include <stdlib.h> @@ -51,7 +75,129 @@ #include "BKE_softbody.h" #include "BKE_displist.h" +#include "BIF_editdeform.h" /* ********** soft body engine ******* */ +#define SOFTGOALSNAP 0.999f +// if bp-> goal is above make it a *forced follow original* and skip all ODE stuff for this bp +// removes *unnecessary* stiffnes from ODE system +#define HEUNWARNLIMIT 1 // 50 would be fine i think for detecting severe *stiff* stuff + +float SoftHeunTol = 1.0f; // humm .. this should be calculated from sb parameters and sizes +float steptime = 1.0f/25.0f; // translate framerate to *real* time +float rescale_grav_to_framerate = 1.0f; // since unit of g is [m/sec^2] we need translation from frames to physics time +float rescale_friction_to_framerate = 1.0f; // since unit of drag is [kg/sec] we need translation from frames to physics time + +short SB_ENABLE = 0; // quick hack to switch sb integration in 3d header + + +void softbody_scale_time(float steptime) +{ + rescale_grav_to_framerate = steptime*steptime; + rescale_friction_to_framerate = steptime; +} + + +static int count_quads( Mesh *me) +{ + int a,result = 0; + MFace *mface= me->mface; + if(mface ) { + for(a=me->totface; a>0; a--, mface++) {if(mface->v4) result++;} + } + return result; +} + +static void add_quad_diag_springs(Object *ob) +{ + Mesh *me= ob->data; + MFace *mface= me->mface; + BodyPoint *bp; + BodySpring *bs, *bs_new; + int a ; + + if (ob->soft){ + int nofquads; + nofquads = count_quads(me); + if (nofquads) { + /* resize spring-array to hold additional quad springs */ + bs_new= MEM_callocN( (ob->soft->totspring + nofquads *2 )*sizeof(BodySpring), "bodyspring"); + memcpy(bs_new,ob->soft->bspring,(ob->soft->totspring )*sizeof(BodySpring)); + MEM_freeN(ob->soft->bspring); /* do this before reassigning the pointer or have a 1st class memory leak */ + ob->soft->bspring = bs_new; + /* fill the tail */ + a = 0; + bs = bs_new+ob->soft->totspring; + bp= ob->soft->bpoint; + if(mface ) { + for(a=me->totface; a>0; a--, mface++) { + if(mface->v4) + { + bs->v1= mface->v1; + bs->v2= mface->v3; + bs->strength= 1.0; + bs->len= VecLenf( (bp+bs->v1)->origS, (bp+bs->v2)->origS); + bs++; + bs->v1= mface->v2; + bs->v2= mface->v4; + bs->strength= 1.0; + bs->len= VecLenf( (bp+bs->v1)->origS, (bp+bs->v2)->origS); + bs++; + + } + } + } + + /* now we can announce new springs */ + ob->soft->totspring += nofquads *2; + + } + } +} + + +static void add_bp_springlist(BodyPoint *bp,int springID) +{ + int *newlist; + if (bp->springs == NULL) { + bp->springs = MEM_callocN( sizeof(int), "bpsprings"); + bp->springs[0] = springID; + bp->nofsprings = 1; + } + else { + bp->nofsprings++; + newlist = MEM_callocN(bp->nofsprings * sizeof(int), "bpsprings"); + memcpy(newlist,bp->springs,(bp->nofsprings-1)* sizeof(int)); + MEM_freeN(bp->springs); + bp->springs = newlist; + bp->springs[bp->nofsprings-1] = springID; + } +} + +/* do this once when sb is build +it is O(N^2) so scanning for springs every iteration is too expensive +*/ +static void build_bps_springlist(Object *ob) +{ + SoftBody *sb= ob->soft; // is supposed to be there + BodyPoint *bp; + BodySpring *bs; + + int a,b; + if (!sb) return; // paranoya check + for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { + /* scan for attached inner springs */ + for(b=sb->totspring, bs= sb->bspring; b>0; b--, bs++) { + if (( (sb->totpoint-a) == bs->v1) ){ + add_bp_springlist(bp,sb->totspring -b); + } + if (( (sb->totpoint-a) == bs->v2) ){ + add_bp_springlist(bp,sb->totspring -b); + } + }//for springs + if (bp->nofsprings) printf(" node %d has %d spring links\n",a,bp->nofsprings); + }//for bp +} + static SoftBody *new_softbody(int totpoint, int totspring) @@ -73,80 +219,361 @@ static SoftBody *new_softbody(int totpoint, int totspring) void free_softbody(SoftBody *sb) { if(sb) { - if(sb->bpoint) MEM_freeN(sb->bpoint); + int a; + BodyPoint *bp; + if(sb->bpoint){ + for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { + /* free spring list */ + if (bp->springs != NULL) { + MEM_freeN(bp->springs); + } + } + MEM_freeN(sb->bpoint); + } if(sb->bspring) MEM_freeN(sb->bspring); - MEM_freeN(sb); } } +/* ************ dynamics ********** */ + +/* aye this belongs to arith.c */ +void Vec3PlusStVec(float *v, float s, float *v1) +{ + v[0] += s*v1[0]; + v[1] += s*v1[1]; + v[2] += s*v1[2]; +} + -static void softbody_calc_forces(Object *ob) +static void softbody_calc_forces(Object *ob, float dtime) { SoftBody *sb= ob->soft; // is supposed to be there - BodyPoint *bp; - float ks; - int a; - + BodyPoint *bp; + float iks,ks,kd,gravity,actspringlen,forcefactor,sd[3]; + int a,b; + BodyPoint *bproot; + BodySpring *bs; /* clear forces */ for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { bp->force[0]= bp->force[1]= bp->force[2]= 0.0; } - - /* spring constant */ - ks= ob->springf; - - /* accumulate forces, vertex stiffness */ + + gravity = ob->sb_nodemass * ob->sb_grav * rescale_grav_to_framerate; + iks = 1.0f/(1.0f-ob->sb_inspring)-1.0f ;/* inner spring constants function */ + bproot= sb->bpoint; /* need this for proper spring addressing */ + for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { - - bp->force[0]= ks*(bp->orig[0]-bp->pos[0]); - bp->force[1]= ks*(bp->orig[1]-bp->pos[1]); - bp->force[2]= ks*(bp->orig[2]-bp->pos[2]); + if(bp->goal < SOFTGOALSNAP){ // ommit this bp when i snaps + float auxvect[3]; // aux unit vector + float velgoal[3]; + float absvel =0, projvel= 0; + + /* do goal stuff */ + /* true elastic goal */ + VecSubf(auxvect,bp->origT,bp->pos); + ks = 1.0f/(1.0f- bp->goal*ob->sb_goalspring)-1.0f ; + bp->force[0]= ks*(auxvect[0]); + bp->force[1]= ks*(auxvect[1]); + bp->force[2]= ks*(auxvect[2]); + /* calulate damping forces generated by goals*/ + VecSubf(velgoal,bp->origS, bp->origE); + kd = ob->sb_goalfrict * rescale_friction_to_framerate ; + + if (dtime > 0.0 ) { // make sure friction does not become rocket motor on time reversal + bp->force[0]-= kd * (velgoal[0] + bp->vec[0]); + bp->force[1]-= kd * (velgoal[1] + bp->vec[1]); + bp->force[2]-= kd * (velgoal[2] + bp->vec[2]); + } + else { + bp->force[0]-= kd * (velgoal[0] - bp->vec[0]); + bp->force[1]-= kd * (velgoal[1] - bp->vec[1]); + bp->force[2]-= kd * (velgoal[2] - bp->vec[2]); + } + /* done goal stuff */ + + + /* gravitation */ + bp->force[2]-= gravity*ob->sb_nodemass; /* individual mass of node here */ + + /* friction in media */ + kd= ob->sb_mediafrict* rescale_friction_to_framerate; + /* assume it to be proportional to actual velocity */ + bp->force[0]-= bp->vec[0]*kd; + bp->force[1]-= bp->vec[1]*kd; + bp->force[2]-= bp->vec[2]*kd; + /* friction in media done */ + + /*other forces*/ + /* this is the place where other forces can be added + yes, constraints and collision stuff should go here too (read baraff papers on that!) + */ + /*other forces done*/ + + /* nice things could be done with anisotropic friction + like wind/air resistance in normal direction + --> having a piece of cloth sailing down + but this needs to have a *valid* vertex normal + *valid* means to be calulated on time axis + hrms .. may be a rough one could be used as well .. let's see + */ + + if (1){ /* big mesh optimization */ + /* run over attached inner spring list */ + if (sb->bspring){ // spring list exists at all ? + for(b=bp->nofsprings;b>0;b--){ + bs = sb->bspring + bp->springs[b-1]; + if (( (sb->totpoint-a) == bs->v1) ){ + actspringlen= VecLenf( (bproot+bs->v2)->pos, bp->pos); + VecSubf(sd,(bproot+bs->v2)->pos, bp->pos); + Normalise(sd); + + // friction stuff V1 + VecSubf(velgoal,bp->vec,(bproot+bs->v2)->vec); + kd = ob->sb_infrict * rescale_friction_to_framerate ; + absvel = Normalise(velgoal); + projvel = ABS(Inpf(sd,velgoal)); + kd *= absvel * projvel; + Vec3PlusStVec(bp->force,-kd,velgoal); + + if(bs->len > 0.0) /* check for degenerated springs */ + forcefactor = (bs->len - actspringlen)/bs->len * iks; + else + forcefactor = actspringlen * iks; + + Vec3PlusStVec(bp->force,-forcefactor,sd); + + } + + if (( (sb->totpoint-a) == bs->v2) ){ + actspringlen= VecLenf( (bproot+bs->v1)->pos, bp->pos); + VecSubf(sd,bp->pos,(bproot+bs->v1)->pos); + Normalise(sd); + + // friction stuff V2 + VecSubf(velgoal,bp->vec,(bproot+bs->v1)->vec); + kd = ob->sb_infrict * rescale_friction_to_framerate ; + absvel = Normalise(velgoal); + projvel = ABS(Inpf(sd,velgoal)); + kd *= absvel * projvel; + Vec3PlusStVec(bp->force,-kd,velgoal); + + if(bs->len > 0.0) + forcefactor = (bs->len - actspringlen)/bs->len * iks; + else + forcefactor = actspringlen * iks; + Vec3PlusStVec(bp->force,+forcefactor,sd); + } + } + } //if spring list exists at all ? + } + else{ // this branch is not completly uptaded for friction stuff + /* scan for attached inner springs makes it a O(N^2) thing = bad !*/ + /* obsolete .. but if someone wants to try the effect :) */ + for(b=sb->totspring, bs= sb->bspring; b>0; b--, bs++) { + if (( (sb->totpoint-a) == bs->v1) ){ + actspringlen= VecLenf( (bproot+bs->v2)->pos, bp->pos); + VecSubf(sd,(bproot+bs->v2)->pos, bp->pos); + Normalise(sd); + + + if(bs->len > 0.0) /* check for degenerated springs */ + forcefactor = (bs->len - actspringlen)/bs->len * iks; + else + forcefactor = actspringlen * iks; + Vec3PlusStVec(bp->force,-forcefactor,sd); + } + + if (( (sb->totpoint-a) == bs->v2) ){ + actspringlen= VecLenf( (bproot+bs->v1)->pos, bp->pos); + VecSubf(sd,bp->pos,(bproot+bs->v1)->pos); + Normalise(sd); + + if(bs->len > 0.0) + forcefactor = (bs->len - actspringlen)/bs->len * iks; + else + forcefactor = actspringlen * iks; + Vec3PlusStVec(bp->force,+forcefactor,sd); + } + }// no snap + }//for + } } } -static void softbody_apply_forces(Object *ob, float dtime) +static void softbody_apply_forces(Object *ob, float dtime, int mode, float *err) { + /* time evolution */ + /* actually does an explicit euler step mode == 0 */ + /* or heun ~ 2nd order runge-kutta steps, mode 1,2 */ SoftBody *sb= ob->soft; // is supposed to be there BodyPoint *bp; - float kd; + float dx[3],dv[3]; int a; + float timeovermass; + float maxerr = 0.0; - kd= 1.0-ob->damping; + if (ob->sb_nodemass > 0.09999f) timeovermass = dtime/ob->sb_nodemass; + else timeovermass = dtime/0.09999f; for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { - // friction - bp->vec[0]*= kd; - bp->vec[1]*= kd; - bp->vec[2]*= kd; - - VECADD(bp->vec, bp->vec, bp->force); // mass here!? - VECADD(bp->pos, bp->pos, bp->vec); + if(bp->goal < SOFTGOALSNAP){ + + /* so here is dv/dt = a = sum(F_springs)/m + gravitation + some friction forces */ + /* the euler step for velocity then becomes */ + /* v(t + dt) = v(t) + a(t) * dt */ + bp->force[0]*= timeovermass; /* individual mass of node here */ + bp->force[1]*= timeovermass; + bp->force[2]*= timeovermass; + /* some nasty if's to have heun in here too */ + VECCOPY(dv,bp->force); + if (mode == 1){ + VECCOPY(bp->prevvec,bp->vec); + VECCOPY(bp->prevdv ,dv); + } + if (mode ==2){ + /* be optimistic and execute step */ + bp->vec[0] = bp->prevvec[0] + 0.5f * (dv[0] + bp->prevdv[0]); + bp->vec[1] = bp->prevvec[1] + 0.5f * (dv[1] + bp->prevdv[1]); + bp->vec[2] = bp->prevvec[2] + 0.5f * (dv[2] + bp->prevdv[2]); + /* compare euler to heun to estimate error for step sizing */ + maxerr = MAX2(maxerr,ABS(dv[0] - bp->prevdv[0])); + maxerr = MAX2(maxerr,ABS(dv[1] - bp->prevdv[1])); + maxerr = MAX2(maxerr,ABS(dv[2] - bp->prevdv[2])); + } + else {VECADD(bp->vec, bp->vec, bp->force);} + + /* so here is dx/dt = v */ + /* the euler step for location then becomes */ + /* x(t + dt) = x(t) + v(t) * dt */ + + VECCOPY(dx,bp->vec); + dx[0]*=dtime ; + dx[1]*=dtime ; + dx[2]*=dtime ; + + /* again some nasty if's to have heun in here too */ + if (mode ==1){ + VECCOPY(bp->prevpos,bp->pos); + VECCOPY(bp->prevdx ,dx); + } + + if (mode ==2){ + bp->pos[0] = bp->prevpos[0] + 0.5f * ( dx[0] + bp->prevdx[0]); + bp->pos[1] = bp->prevpos[1] + 0.5f * ( dx[1] + bp->prevdx[1]); + bp->pos[2] = bp->prevpos[2] + 0.5f* ( dx[2] + bp->prevdx[2]); + maxerr = MAX2(maxerr,ABS(dx[0] - bp->prevdx[0])); + maxerr = MAX2(maxerr,ABS(dx[1] - bp->prevdx[1])); + maxerr = MAX2(maxerr,ABS(dx[2] - bp->prevdx[2])); + } + else { VECADD(bp->pos, bp->pos, dx);} + + }//snap + } //for + if (err){ /* so step size will be controlled by biggest difference in slope */ + *err = maxerr; } } +/* used by heun when it overshoots */ +static void softbody_restore_prev_step(Object *ob) +{ + SoftBody *sb= ob->soft; // is supposed to be there + BodyPoint *bp; + int a; + for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { + VECCOPY(bp->vec,bp->prevvec); + VECCOPY(bp->pos,bp->prevpos); + } +} + + +/* unused */ static void softbody_apply_goal(Object *ob, float dtime) { + SoftBody *sb= ob->soft; // is supposed to be there BodyPoint *bp; float vec[3], ks; int a; for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { - ks= bp->goal; + ks= bp->goal*dtime; // this is hackish, screws up physics but stabilizes - vec[0]= ks*(bp->orig[0]-bp->pos[0]); - vec[1]= ks*(bp->orig[1]-bp->pos[1]); - vec[2]= ks*(bp->orig[2]-bp->pos[2]); + vec[0]= ks*(bp->origT[0]-bp->pos[0]); + vec[1]= ks*(bp->origT[1]-bp->pos[1]); + vec[2]= ks*(bp->origT[2]-bp->pos[2]); VECADD(bp->pos, bp->pos, vec); - ks= 1.0-ks; + ks= 1.0f-ks; bp->vec[0]*= ks; bp->vec[1]*= ks; bp->vec[2]*= ks; + + } +} + +static void softbody_apply_goalsnap(Object *ob) +{ + SoftBody *sb= ob->soft; // is supposed to be there + BodyPoint *bp; + int a; + + for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { + if (bp->goal >= SOFTGOALSNAP){ + VECCOPY(bp->prevpos,bp->pos); + VECCOPY(bp->pos,bp->origT); + } + } +} + +static void softbody_force_goal(Object *ob) +{ + SoftBody *sb= ob->soft; // is supposed to be there + BodyPoint *bp; + int a; + + for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) { + VECCOPY(bp->pos,bp->origT); + bp->vec[0] = bp->origE[0] - bp->origS[0]; + bp->vec[1] = bp->origE[1] - bp->origS[1]; + bp->vec[2] = bp->origE[2] - bp->origS[2]; } } +static void interpolate_exciter(Object *ob, int timescale, int time) +{ + Mesh *me= ob->data; + //MEdge *medge= me->medge; + int a; + BodyPoint *bp; + float f; + + if(ob->soft) { + f = (float)time/(float)timescale; + bp= ob->soft->bpoint; + for(a=0; a<me->totvert; a++, bp++) { + bp->origT[0] = bp->origS[0] + f*(bp->origE[0] - bp->origS[0]); + bp->origT[1] = bp->origS[1] + f*(bp->origE[1] - bp->origS[1]); + bp->origT[2] = bp->origS[2] + f*(bp->origE[2] - bp->origS[2]); + if (bp->goal >= SOFTGOALSNAP){ + bp->vec[0] = bp->origE[0] - bp->origS[0]; + bp->vec[1] = bp->origE[1] - bp->origS[1]; + bp->vec[2] = bp->origE[2] - bp->origS[2]; + } + } + /* hrms .. do springs alter their lenght ? + if(medge) { + bs= ob->soft->bspring; + bp= ob->soft->bpoint; + for(a=0; (a<me->totedge && a < ob->soft->totspring ); a++, medge++, bs++) { + bs->len= VecLenf( (bp+bs->v1)->origT, (bp+bs->v2)->origT); + } + } + */ + } +} + /* ************ convertors ********** */ @@ -157,27 +584,60 @@ static void mesh_update_softbody(Object *ob) MVert *mvert= me->mvert; MEdge *medge= me->medge; BodyPoint *bp; - BodySpring *bs; int a; - if(ob->soft) { bp= ob->soft->bpoint; for(a=0; a<me->totvert; a++, mvert++, bp++) { - VECCOPY(bp->orig, mvert->co); - Mat4MulVecfl(ob->obmat, bp->orig); + VECCOPY(bp->origS, bp->origE); + VECCOPY(bp->origE, mvert->co); + Mat4MulVecfl(ob->obmat, bp->origE); + VECCOPY(bp->origT, bp->origE); } + /* hrms .. do springs alter their lenght ? if(medge) { bs= ob->soft->bspring; bp= ob->soft->bpoint; - for(a=0; a<me->totedge; a++, medge++, bs++) { - bs->len= VecLenf( (bp+bs->v1)->orig, (bp+bs->v2)->orig); + for(a=0; (a<me->totedge && a < ob->soft->totspring ); a++, medge++, bs++) { + bs->len= VecLenf( (bp+bs->v1)->origE, (bp+bs->v2)->origE); } } + */ } } + +int get_scalar_from_named_vertexgroup(Object *ob, char *name, int vertID, float *target) +/* result 0 on success, else indicates error number +-- kind of *inverse* result defintion, +-- but this way we can signal error condition to caller +-- and yes this function must not be here but in a *vertex group module* +*/ +{ + int i,groupindex; + bDeformGroup *locGroup = NULL; + MDeformVert *dv; + locGroup=get_named_vertexgroup (ob,name); + if(locGroup){ + /* retrieve index for that group */ + groupindex = get_defgroup_num (ob,locGroup); + /* spot the vert in deform vert list at mesh */ + /* todo (coder paranoya) what if ob->data is not a mesh .. */ + /* hrms.. would like to have the same for lattices anyhoo */ + dv = ((Mesh*)ob->data)->dvert + vertID; + /* Lets see if this vert is in the weight group */ + for (i=0; i<dv->totweight; i++){ + if (dv->dw[i].def_nr == groupindex){ + *target=dv->dw[i].weight; /* got it ! */ + return 0; + } + }/*for*/ + return 2; + }/*if(locGroup)*/ + return 1; +} + /* makes totally fresh start situation */ static void mesh_to_softbody(Object *ob) { @@ -195,11 +655,44 @@ static void mesh_to_softbody(Object *ob) for(a=me->totvert; a>0; a--, mvert++, bp++) { VECCOPY(bp->pos, mvert->co); Mat4MulVecfl(ob->obmat, bp->pos); // yep, sofbody is global coords - VECCOPY(bp->orig, bp->pos); + VECCOPY(bp->origS, bp->pos); + VECCOPY(bp->origE, bp->pos); + VECCOPY(bp->origT, bp->pos); bp->vec[0]= bp->vec[1]= bp->vec[2]= 0.0; bp->weight= 1.0; bp->goal= 0.5; + bp->nofsprings=0; + bp->springs=NULL; + if (1) { /* switch to vg scalars*/ + /* get scalar values needed *per vertex* from vertex group functions, + so we can *paint* them nicly .. + they are normalized [0.0..1.0] so may be we need amplitude for scale + which can be done by caller + but still .. i'd like it to go this way + */ + int error; + char name[32] = "SOFTGOAL"; + float temp; + error = get_scalar_from_named_vertexgroup(ob,name,me->totvert - a,&temp); + if (!error) bp->goal = temp; + if (bp->goal < ob->sb_mingoal) bp->goal = ob->sb_mingoal; + if (bp->goal > ob->sb_maxgoal) bp->goal = ob->sb_maxgoal; + /* a little ad hoc changing the goal control to be less *sharp* */ + bp->goal = (float)pow(bp->goal,4.0f); +/* to proove the concept +this would enable per vertex *mass painting* + strcpy(name,"SOFTMASS"); + error = get_scalar_from_named_vertexgroup(ob,name,me->totvert - a,&temp); + if (!error) bp->mass = temp * ob->rangeofmass; +*/ + + + + } /* switch to vg scalars */ } + + + if(medge) { bs= ob->soft->bspring; bp= ob->soft->bpoint; @@ -207,10 +700,20 @@ static void mesh_to_softbody(Object *ob) bs->v1= medge->v1; bs->v2= medge->v2; bs->strength= 1.0; - bs->len= VecLenf( (bp+bs->v1)->orig, (bp+bs->v2)->orig); + bs->len= VecLenf( (bp+bs->v1)->origS, (bp+bs->v2)->origS); } } + + /* insert *diagonal* springs in quads if desired */ + if (ob->softflag & 0x02) { + add_quad_diag_springs(ob); + + } + + build_bps_springlist(ob); /* big mesh optimization */ + /* vertex colors are abused as weights here, however they're stored in faces... uhh */ + /* naah .. we don't do it any more bjornmose :-) if(mface && me->mcol) { char *mcol= (char *)me->mcol; for(a=me->totface; a>0; a--, mface++, mcol+=16) { @@ -234,6 +737,7 @@ static void mesh_to_softbody(Object *ob) } } } + */ bp= ob->soft->bpoint; for(a=me->totvert; a>0; a--, bp++) { //printf("a %d goal %f\n", a, bp->goal); @@ -292,7 +796,7 @@ void object_update_softbody(Object *ob) } /* makes totally fresh start situation */ -void object_to_softbody(Object *ob) +void object_to_softbody(Object *ob,float ctime) { if(ob->soft) free_softbody(ob->soft); @@ -301,6 +805,7 @@ void object_to_softbody(Object *ob) switch(ob->type) { case OB_MESH: mesh_to_softbody(ob); + ob->soft->ctime = ctime; break; case OB_LATTICE: lattice_to_softbody(ob); @@ -326,32 +831,165 @@ void softbody_to_object(Object *ob) /* simulates one step. ctime is in frames not seconds */ + void object_softbody_step(Object *ob, float ctime) { float dtime; - + int timescale,t; + float forcetime; + float err; + + /* this is a NO! NO! + ========================== if(ob->soft==NULL) { object_to_softbody(ob); if(ob->soft==NULL) return; ob->soft->ctime= ctime; } + // you can't create a soft object on the fly + // 1. inner spings need a *default* length for crinkles/wrinkles, + // surface area and volume preservation + // 2. initial conditions for velocities and positions need to be defined + // for a certain point of time .. say t0 + // 3. and since we have friction and *outer* movement + // the history of the *outer* movements will affect where we end up + // sooo atm going to edit mode and back ( back calls object_to_softbody(ob,1.0f) + is the only way to create softbody data + */ + + /* first attempt to set initial conditions for softbodies + rules + 1. ODE solving is disabled / via button in 3dview header /otherways do regular softbody stuff + 2. set SB positions to *goal* + 3. set SB velocities to match *goal* movement + + */ + if (SB_ENABLE == 0){ + if(ob->soft==NULL) { + return; /* nothing to do */ + } + object_update_softbody(ob); + ob->soft->ctime= ctime; + interpolate_exciter(ob,200,200); + softbody_force_goal(ob); + softbody_to_object(ob); + return; /* no dynamics wanted */ + + } + + if(ob->soft==NULL) { + /* aye no soft object created bail out here */ + printf("Softbody Zombie \n"); + return; + } + softbody_scale_time(steptime); // translate frames/sec and lenghts unit to SI system dtime= ctime - ob->soft->ctime; - dtime= ABS(dtime); - if(dtime > 0.0) { - /* desired vertex locations in oldloc */ + // dtime= ABS(dtime); no no we want to go back in time with IPOs + timescale = (int)(ob->softtime * ABS(dtime)); + if(ABS(dtime) > 0.0) { object_update_softbody(ob); + if (ob->softflag & 0x04){ + /* special case of 2nd order Runge-Kutta type AKA Heun */ + float timedone =0.0; + /* counter for emergency brake + * we don't want to lock up the system if physics fail + */ + int loops =0 ; + SoftHeunTol = ob->softtime; // humm .. this should be calculated from sb parameters and sizes + + forcetime = dtime; /* hope for integrating in one step */ + while ( (ABS(timedone) < ABS(dtime)) && (loops < 2000) ) + { + if (ABS(dtime) > 3.0 ){ + printf("SB_STEPSIZE \n"); + break; // sorry but i must assume goal movement can't be interpolated any more + } + //set goals in time + interpolate_exciter(ob,200,(int)(200.0*(timedone/dtime))); + // do predictive euler step + softbody_calc_forces(ob,forcetime); + softbody_apply_forces(ob,forcetime,1, NULL); + // crop new slope values to do averaged slope step + softbody_calc_forces(ob,forcetime); + softbody_apply_forces(ob,forcetime,2, &err); + softbody_apply_goalsnap(ob); + + if (err > SoftHeunTol){ // error needs to be scaled to some quantity + softbody_restore_prev_step(ob); + forcetime /= 2.0; + } + else { + + float newtime = forcetime * 1.1f; // hope for 1.1 times better conditions in next step + if (err > SoftHeunTol/2.0){ // stay with this stepsize unless err really small + newtime = forcetime; + } + timedone += forcetime; + if (forcetime > 0.0) + forcetime = MIN2(dtime - timedone,newtime); + else + forcetime = MAX2(dtime - timedone,newtime); + } + loops++; + } + // move snapped to final position + interpolate_exciter(ob,2,2); + softbody_apply_goalsnap(ob); + if (loops > HEUNWARNLIMIT) /* monitor high loop counts say 1000 after testing */ + printf("%d heun integration loops/frame \n",loops); + } + else + /* do brute force explicit euler */ + /* inner intagration loop */ + /* */ + // loop n times so that n*h = duration of one frame := 1 + // x(t+h) = x(t) + h*v(t); + // v(t+h) = v(t) + h*f(x(t),t); + for(t=1 ; t <= timescale; t++) { + if (ABS(dtime) > 15 ) break; + + /* the *goal* mesh must use the n*h timing too ! + use *cheap* linear intepolation for that */ + interpolate_exciter(ob,timescale,t); + if (timescale > 0 ) + { + forcetime = dtime/timescale; + + /* does not fit the concept sloving ODEs :) */ + /* softbody_apply_goal(ob,forcetime ); */ + + /* explicit Euler integration */ + /* we are not controling a nuclear power plant! + so rought *almost* physical behaviour is acceptable. + in cases of *mild* stiffnes cranking up timscale -> decreasing stepsize *h* + avoids instability */ + softbody_calc_forces(ob,forcetime); + softbody_apply_forces(ob,forcetime,0, NULL); + softbody_apply_goalsnap(ob); - /* extra for desired vertex locations */ - softbody_apply_goal(ob, dtime); - - softbody_calc_forces(ob); - softbody_apply_forces(ob, dtime); +// if (0){ + /* ok here comes the überhammer + use a semi implicit euler integration to tackle *all* stiff conditions + but i doubt the cost/benifit holds for most of the cases + -- to be coded*/ +// } + + } + } + /* and apply to vertices */ - softbody_to_object(ob); + softbody_to_object(ob); ob->soft->ctime= ctime; + } // if(ABS(dtime) > 0.0) + else { + // rule : you have asked for the current state of the softobject + // since dtime= ctime - ob->soft->ctime; + // and we were not notifified about any other time changes + // so here it is ! + softbody_to_object(ob); } } |