5
* ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
7
* This program is free software; you can redistribute it and/or
8
* modify it under the terms of the GNU General Public License
9
* as published by the Free Software Foundation; either version 2
10
* of the License, or (at your option) any later version. The Blender
11
* Foundation also sells licenses for use in proprietary software under
12
* the Blender License. See http://www.blender.org/BL/ for information
15
* This program is distributed in the hope that it will be useful,
16
* but WITHOUT ANY WARRANTY; without even the implied warranty of
17
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18
* GNU General Public License for more details.
20
* You should have received a copy of the GNU General Public License
21
* along with this program; if not, write to the Free Software Foundation,
22
* Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
24
* The Original Code is Copyright (C) Blender Foundation
25
* All rights reserved.
27
* The Original Code is: all of this file.
29
* Contributor(s): none yet.
31
* ***** END GPL/BL DUAL LICENSE BLOCK *****
36
variables on the UI for now
38
float mediafrict; friction to env
39
float nodemass; softbody mass of *vertex*
40
float grav; softbody amount of gravitaion to apply
42
float goalspring; softbody goal springs
43
float goalfrict; softbody goal springs friction
44
float mingoal; quick limits for goal
47
float inspring; softbody inner springs
48
float infrict; softbody inner springs friction
58
#include "MEM_guardedalloc.h"
61
#include "DNA_curve_types.h"
62
#include "DNA_object_types.h"
63
#include "DNA_object_force.h" /* here is the softbody struct */
64
#include "DNA_key_types.h"
65
#include "DNA_mesh_types.h"
66
#include "DNA_meshdata_types.h"
67
#include "DNA_lattice_types.h"
68
#include "DNA_scene_types.h"
70
#include "BLI_blenlib.h"
71
#include "BLI_arithb.h"
73
#include "BKE_displist.h"
74
#include "BKE_effect.h"
75
#include "BKE_global.h"
77
#include "BKE_object.h"
78
#include "BKE_softbody.h"
79
#include "BKE_utildefines.h"
81
#include "BIF_editdeform.h"
83
extern bDeformGroup *get_named_vertexgroup(Object *ob, char *name);
84
extern int get_defgroup_num (Object *ob, bDeformGroup *dg);
87
/* ********** soft body engine ******* */
88
#define SOFTGOALSNAP 0.999f
89
// if bp-> goal is above make it a *forced follow original* and skip all ODE stuff for this bp
90
// removes *unnecessary* stiffnes from ODE system
91
#define HEUNWARNLIMIT 1 // 50 would be fine i think for detecting severe *stiff* stuff
94
float SoftHeunTol = 1.0f; // humm .. this should be calculated from sb parameters and sizes
96
/* local prototypes */
97
static void free_softbody_intern(SoftBody *sb);
100
/*+++ frame based timing +++*/
102
//physical unit of force is [kg * m / sec^2]
104
static float sb_grav_force_scale(Object *ob)
105
// since unit of g is [m/sec^2] and F = mass * g we rescale unit mass of node to 1 gramm
106
// put it to a function here, so we can add user options later without touching simulation code
111
static float sb_fric_force_scale(Object *ob)
112
// rescaling unit of drag [1 / sec] to somehow reasonable
113
// put it to a function here, so we can add user options later without touching simulation code
118
static float sb_time_scale(Object *ob)
119
// defining the frames to *real* time relation
121
SoftBody *sb= ob->soft; // is supposed to be there
123
return(sb->physics_speed); //hrms .. this could be IPO as well :)
124
// estimated range [0.001 sluggish slug - 100.0 very fast (i hope ODE solver can handle that)]
125
// 1 approx = a unit 1 pendulum at g = 9.8 [earth conditions] has period 65 frames
126
// theory would give a 50 frames period .. so there must be something inaccurate .. looking for that (BM)
130
this would be frames/sec independant timing assuming 25 fps is default
131
but does not work very well with NLA
132
return (25.0f/G.scene->r.frs_sec)
135
/*--- frame based timing ---*/
138
static int count_mesh_quads(Mesh *me)
141
MFace *mface= me->mface;
144
for(a=me->totface; a>0; a--, mface++) {
145
if(mface->v4) result++;
151
static void add_mesh_quad_diag_springs(Object *ob)
154
MFace *mface= me->mface;
156
BodySpring *bs, *bs_new;
162
nofquads = count_mesh_quads(me);
164
/* resize spring-array to hold additional quad springs */
165
bs_new= MEM_callocN( (ob->soft->totspring + nofquads *2 )*sizeof(BodySpring), "bodyspring");
166
memcpy(bs_new,ob->soft->bspring,(ob->soft->totspring )*sizeof(BodySpring));
168
if(ob->soft->bspring)
169
MEM_freeN(ob->soft->bspring); /* do this before reassigning the pointer or have a 1st class memory leak */
170
ob->soft->bspring = bs_new;
174
bs = bs_new+ob->soft->totspring;
175
bp= ob->soft->bpoint;
177
for(a=me->totface; a>0; a--, mface++) {
182
bs->len= VecLenf( (bp+bs->v1)->origS, (bp+bs->v2)->origS);
187
bs->len= VecLenf( (bp+bs->v1)->origS, (bp+bs->v2)->origS);
194
/* now we can announce new springs */
195
ob->soft->totspring += nofquads *2;
201
static void add_bp_springlist(BodyPoint *bp,int springID)
205
if (bp->springs == NULL) {
206
bp->springs = MEM_callocN( sizeof(int), "bpsprings");
207
bp->springs[0] = springID;
212
newlist = MEM_callocN(bp->nofsprings * sizeof(int), "bpsprings");
213
memcpy(newlist,bp->springs,(bp->nofsprings-1)* sizeof(int));
214
MEM_freeN(bp->springs);
215
bp->springs = newlist;
216
bp->springs[bp->nofsprings-1] = springID;
220
/* do this once when sb is build
221
it is O(N^2) so scanning for springs every iteration is too expensive
223
static void build_bps_springlist(Object *ob)
225
SoftBody *sb= ob->soft; // is supposed to be there
230
if (sb==NULL) return; // paranoya check
232
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
233
/* scan for attached inner springs */
234
for(b=sb->totspring, bs= sb->bspring; b>0; b--, bs++) {
235
if (( (sb->totpoint-a) == bs->v1) ){
236
add_bp_springlist(bp,sb->totspring -b);
238
if (( (sb->totpoint-a) == bs->v2) ){
239
add_bp_springlist(bp,sb->totspring -b);
242
// if (bp->nofsprings) printf(" node %d has %d spring links\n",a,bp->nofsprings);
247
/* creates new softbody if didn't exist yet, makes new points and springs arrays */
248
/* called in mesh_to_softbody */
249
static void renew_softbody(Object *ob, int totpoint, int totspring)
253
if(ob->soft==NULL) ob->soft= sbNew();
254
else free_softbody_intern(ob->soft);
258
sb->totpoint= totpoint;
259
sb->totspring= totspring;
261
sb->bpoint= MEM_mallocN( totpoint*sizeof(BodyPoint), "bodypoint");
263
sb->bspring= MEM_mallocN( totspring*sizeof(BodySpring), "bodyspring");
267
static void free_softbody_baked(SoftBody *sb)
272
for(k=0; k<sb->totkey; k++) {
273
key= *(sb->keys + k);
274
if(key) MEM_freeN(key);
276
if(sb->keys) MEM_freeN(sb->keys);
283
/* only frees internal data */
284
static void free_softbody_intern(SoftBody *sb)
291
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
292
/* free spring list */
293
if (bp->springs != NULL) {
294
MEM_freeN(bp->springs);
297
MEM_freeN(sb->bpoint);
300
if(sb->bspring) MEM_freeN(sb->bspring);
302
sb->totpoint= sb->totspring= 0;
306
free_softbody_baked(sb);
311
/* ************ dynamics ********** */
313
/* aye this belongs to arith.c */
314
static void Vec3PlusStVec(float *v, float s, float *v1)
321
static int sb_deflect_face(Object *ob,float *actpos, float *futurepos,float *collisionpos, float *facenormal,float *force,float *cf ,float *bounce)
324
float s_actpos[3], s_futurepos[3];
325
VECCOPY(s_actpos,actpos);
327
VECCOPY(s_futurepos,futurepos);
328
if (bounce) *bounce *= 1.5f;
331
deflected= SoftBodyDetectCollision(s_actpos, s_futurepos, collisionpos,
332
facenormal, cf, force , 1,
333
G.scene->r.cfra, ob->lay, ob);
338
/* for future use (BM)
339
static int sb_deflect_edge_face(Object *ob,float *actpos, float *futurepos,float *collisionpos, float *facenormal,float *slip ,float *bounce)
342
float dummy[3],s_actpos[3], s_futurepos[3];
343
SoftBody *sb= ob->soft; // is supposed to be there
344
VECCOPY(s_actpos,actpos);
345
VECCOPY(s_futurepos,futurepos);
346
if (slip) *slip *= 0.98f;
347
if (bounce) *bounce *= 1.5f;
350
deflected= SoftBodyDetectCollision(s_actpos, s_futurepos, collisionpos,
351
facenormal, dummy, dummy , 2,
352
G.scene->r.cfra, ob->lay, ob);
357
// some functions removed here .. to help HOS on next merge (BM)
360
#define USES_DEFLECT 2
361
static int is_there_deflection(unsigned int layer)
366
for(base = G.scene->base.first; base; base= base->next) {
367
if( (base->lay & layer) && base->object->pd) {
368
if(base->object->pd->forcefield) retval |= USES_FIELD;
369
if(base->object->pd->deflect) retval |= USES_DEFLECT;
375
static void softbody_calc_forces(Object *ob, float forcetime)
377
/* rule we never alter free variables :bp->vec bp->pos in here !
378
* this will ruin adaptive stepsize AKA heun! (BM)
380
SoftBody *sb= ob->soft; // is supposed to be there
384
float iks, ks, kd, gravity, actspringlen, forcefactor, sd[3],
385
fieldfactor = 1000.0f,
387
int a, b, do_effector;
390
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
391
bp->force[0]= bp->force[1]= bp->force[2]= 0.0;
394
gravity = sb->grav * sb_grav_force_scale(ob);
396
do_effector= is_there_deflection(ob->lay);
397
iks = 1.0f/(1.0f-sb->inspring)-1.0f ;/* inner spring constants function */
398
bproot= sb->bpoint; /* need this for proper spring addressing */
400
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
401
if(bp->goal < SOFTGOALSNAP){ // ommit this bp when i snaps
402
float auxvect[3]; // aux unit vector
404
float absvel =0, projvel= 0;
407
if(ob->softflag & OB_SB_GOAL) {
408
/* true elastic goal */
409
VecSubf(auxvect,bp->origT,bp->pos);
410
ks = 1.0f/(1.0f- bp->goal*sb->goalspring)-1.0f ;
411
bp->force[0]= ks*(auxvect[0]);
412
bp->force[1]= ks*(auxvect[1]);
413
bp->force[2]= ks*(auxvect[2]);
414
/* calulate damping forces generated by goals*/
415
VecSubf(velgoal,bp->origS, bp->origE);
416
kd = sb->goalfrict * sb_fric_force_scale(ob) ;
418
if (forcetime > 0.0 ) { // make sure friction does not become rocket motor on time reversal
419
bp->force[0]-= kd * (velgoal[0] + bp->vec[0]);
420
bp->force[1]-= kd * (velgoal[1] + bp->vec[1]);
421
bp->force[2]-= kd * (velgoal[2] + bp->vec[2]);
424
bp->force[0]-= kd * (velgoal[0] - bp->vec[0]);
425
bp->force[1]-= kd * (velgoal[1] - bp->vec[1]);
426
bp->force[2]-= kd * (velgoal[2] - bp->vec[2]);
429
/* done goal stuff */
433
bp->force[2]-= gravity*sb->nodemass; /* individual mass of node here */
435
/* particle field & vortex */
436
if(do_effector & USES_FIELD) {
437
float force[3]= {0.0f, 0.0f, 0.0f};
438
float speed[3]= {0.0f, 0.0f, 0.0f};
439
float eval_sb_fric_force_scale = sb_fric_force_scale(ob); // just for calling functio once
441
pdDoEffector(bp->pos, force, speed, (float)G.scene->r.cfra, ob->lay,PE_WIND_AS_SPEED);
442
// note: now we have wind as motion of media, so we can do anisotropic stuff here,
443
// if we had vertex normals here(BM)
444
/* apply forcefield*/
445
VecMulf(force,fieldfactor* eval_sb_fric_force_scale);
446
VECADD(bp->force, bp->force, force);
448
/* friction in moving media */
449
kd= sb->mediafrict* eval_sb_fric_force_scale;
450
bp->force[0] -= kd * (bp->vec[0] + windfactor*speed[0]/eval_sb_fric_force_scale);
451
bp->force[1] -= kd * (bp->vec[1] + windfactor*speed[1]/eval_sb_fric_force_scale);
452
bp->force[2] -= kd * (bp->vec[2] + windfactor*speed[2]/eval_sb_fric_force_scale);
453
/* now we'll have nice centrifugal effect for vortex */
457
/* friction in media (not) moving*/
458
kd= sb->mediafrict* sb_fric_force_scale(ob);
459
/* assume it to be proportional to actual velocity */
460
bp->force[0]-= bp->vec[0]*kd;
461
bp->force[1]-= bp->vec[1]*kd;
462
bp->force[2]-= bp->vec[2]*kd;
463
/* friction in media done */
467
/* this is the place where other forces can be added
468
yes, constraints and collision stuff should go here too (read baraff papers on that!)
470
/* try to match moving collision targets */
471
/* master switch to turn collision off (BM)*/
473
if(do_effector & USES_DEFLECT) {
474
/*sorry for decl. here i'll move 'em up when WIP is done (BM) */
476
float collisionpos[3],facenormal[3];
484
if (sb_deflect_face(ob,bp->pos, bp->pos, collisionpos, facenormal,defforce,&cf,&bounce)){
485
bp->force[0] += defforce[0]*kd;
486
bp->force[1] += defforce[1]*kd;
487
bp->force[2] += defforce[2]*kd;
488
bp->contactfrict = cf;
491
bp->contactfrict = 0.0f;
496
/*other forces done*/
497
/* nice things could be done with anisotropic friction
498
like wind/air resistance in normal direction
499
--> having a piece of cloth sailing down
500
but this needs to have a *valid* vertex normal
501
*valid* means to be calulated on time axis
502
hrms .. may be a rough one could be used as well .. let's see
505
if(ob->softflag & OB_SB_EDGES) {
506
if (sb->bspring){ // spring list exists at all ?
507
for(b=bp->nofsprings;b>0;b--){
508
bs = sb->bspring + bp->springs[b-1];
509
if (( (sb->totpoint-a) == bs->v1) ){
510
actspringlen= VecLenf( (bproot+bs->v2)->pos, bp->pos);
511
VecSubf(sd,(bproot+bs->v2)->pos, bp->pos);
515
VecSubf(velgoal,bp->vec,(bproot+bs->v2)->vec);
516
kd = sb->infrict * sb_fric_force_scale(ob);
517
absvel = Normalise(velgoal);
518
projvel = ABS(Inpf(sd,velgoal));
519
kd *= absvel * projvel;
520
Vec3PlusStVec(bp->force,-kd,velgoal);
522
if(bs->len > 0.0) /* check for degenerated springs */
523
forcefactor = (bs->len - actspringlen)/bs->len * iks;
525
forcefactor = actspringlen * iks;
527
Vec3PlusStVec(bp->force,-forcefactor,sd);
531
if (( (sb->totpoint-a) == bs->v2) ){
532
actspringlen= VecLenf( (bproot+bs->v1)->pos, bp->pos);
533
VecSubf(sd,bp->pos,(bproot+bs->v1)->pos);
537
VecSubf(velgoal,bp->vec,(bproot+bs->v1)->vec);
538
kd = sb->infrict * sb_fric_force_scale(ob);
539
absvel = Normalise(velgoal);
540
projvel = ABS(Inpf(sd,velgoal));
541
kd *= absvel * projvel;
542
Vec3PlusStVec(bp->force,-kd,velgoal);
545
forcefactor = (bs->len - actspringlen)/bs->len * iks;
547
forcefactor = actspringlen * iks;
548
Vec3PlusStVec(bp->force,+forcefactor,sd);
551
}/* existing spring list */
557
static void softbody_apply_forces(Object *ob, float forcetime, int mode, float *err)
560
/* actually does an explicit euler step mode == 0 */
561
/* or heun ~ 2nd order runge-kutta steps, mode 1,2 */
562
SoftBody *sb= ob->soft; // is supposed to be there
569
forcetime *= sb_time_scale(ob);
571
do_effector= is_there_deflection(ob->lay);
573
// claim a minimum mass for vertex
574
if (sb->nodemass > 0.09999f) timeovermass = forcetime/sb->nodemass;
575
else timeovermass = forcetime/0.09999f;
577
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
578
if(bp->goal < SOFTGOALSNAP){
580
/* so here is (v)' = a(cceleration) = sum(F_springs)/m + gravitation + some friction forces + more forces*/
581
/* the ( ... )' operator denotes derivate respective time */
582
/* the euler step for velocity then becomes */
583
/* v(t + dt) = v(t) + a(t) * dt */
584
bp->force[0]*= timeovermass; /* individual mass of node here */
585
bp->force[1]*= timeovermass;
586
bp->force[2]*= timeovermass;
587
/* some nasty if's to have heun in here too */
588
VECCOPY(dv,bp->force);
590
VECCOPY(bp->prevvec, bp->vec);
591
VECCOPY(bp->prevdv, dv);
594
/* be optimistic and execute step */
595
bp->vec[0] = bp->prevvec[0] + 0.5f * (dv[0] + bp->prevdv[0]);
596
bp->vec[1] = bp->prevvec[1] + 0.5f * (dv[1] + bp->prevdv[1]);
597
bp->vec[2] = bp->prevvec[2] + 0.5f * (dv[2] + bp->prevdv[2]);
598
/* compare euler to heun to estimate error for step sizing */
599
maxerr = MAX2(maxerr,ABS(dv[0] - bp->prevdv[0]));
600
maxerr = MAX2(maxerr,ABS(dv[1] - bp->prevdv[1]));
601
maxerr = MAX2(maxerr,ABS(dv[2] - bp->prevdv[2]));
603
else {VECADD(bp->vec, bp->vec, bp->force);}
605
/* so here is (x)'= v(elocity) */
606
/* the euler step for location then becomes */
607
/* x(t + dt) = x(t) + v(t) * dt */
614
/* again some nasty if's to have heun in here too */
616
VECCOPY(bp->prevpos,bp->pos);
617
VECCOPY(bp->prevdx ,dx);
621
bp->pos[0] = bp->prevpos[0] + 0.5f * ( dx[0] + bp->prevdx[0]);
622
bp->pos[1] = bp->prevpos[1] + 0.5f * ( dx[1] + bp->prevdx[1]);
623
bp->pos[2] = bp->prevpos[2] + 0.5f* ( dx[2] + bp->prevdx[2]);
624
maxerr = MAX2(maxerr,ABS(dx[0] - bp->prevdx[0]));
625
maxerr = MAX2(maxerr,ABS(dx[1] - bp->prevdx[1]));
626
maxerr = MAX2(maxerr,ABS(dx[2] - bp->prevdx[2]));
627
/* kind of hack .. while inside collision target .. make movement more *viscous* */
628
if (bp->contactfrict > 0.0f){
629
bp->vec[0] *= (1.0 - bp->contactfrict);
630
bp->vec[1] *= (1.0 - bp->contactfrict);
631
bp->vec[2] *= (1.0 - bp->contactfrict);
634
else { VECADD(bp->pos, bp->pos, dx);}
635
// experimental particle collision suff was here .. just to help HOS on next merge (BM)
638
if (err){ /* so step size will be controlled by biggest difference in slope */
643
/* used by heun when it overshoots */
644
static void softbody_restore_prev_step(Object *ob)
646
SoftBody *sb= ob->soft; // is supposed to be there
650
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
651
VECCOPY(bp->vec, bp->prevvec);
652
VECCOPY(bp->pos, bp->prevpos);
657
static void softbody_apply_goalsnap(Object *ob)
659
SoftBody *sb= ob->soft; // is supposed to be there
663
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
664
if (bp->goal >= SOFTGOALSNAP){
665
VECCOPY(bp->prevpos,bp->pos);
666
VECCOPY(bp->pos,bp->origT);
673
static void softbody_force_goal(Object *ob)
675
SoftBody *sb= ob->soft; // is supposed to be there
679
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
680
VECCOPY(bp->pos,bp->origT);
681
bp->vec[0] = bp->origE[0] - bp->origS[0];
682
bp->vec[1] = bp->origE[1] - bp->origS[1];
683
bp->vec[2] = bp->origE[2] - bp->origS[2];
688
/* expects full initialized softbody */
689
static void interpolate_exciter(Object *ob, int timescale, int time)
691
SoftBody *sb= ob->soft;
696
// note: i removed Mesh usage here, softbody should remain generic! (ton)
698
f = (float)time/(float)timescale;
700
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
701
bp->origT[0] = bp->origS[0] + f*(bp->origE[0] - bp->origS[0]);
702
bp->origT[1] = bp->origS[1] + f*(bp->origE[1] - bp->origS[1]);
703
bp->origT[2] = bp->origS[2] + f*(bp->origE[2] - bp->origS[2]);
704
if (bp->goal >= SOFTGOALSNAP){
705
bp->vec[0] = bp->origE[0] - bp->origS[0];
706
bp->vec[1] = bp->origE[1] - bp->origS[1];
707
bp->vec[2] = bp->origE[2] - bp->origS[2];
711
if(ob->softflag & OB_SB_EDGES) {
712
/* hrms .. do springs alter their lenght ?
713
bs= ob->soft->bspring;
714
bp= ob->soft->bpoint;
715
for(a=0; (a<me->totedge && a < ob->soft->totspring ); a++, bs++) {
716
bs->len= VecLenf( (bp+bs->v1)->origT, (bp+bs->v2)->origT);
723
/* ************ convertors ********** */
725
/* for each object type we need;
726
- xxxx_to_softbody(Object *ob) : a full (new) copy
727
- xxxx_update_softbody(Object *ob) : update refreshes current positions
728
- softbody_to_xxxx(Object *ob) : after simulation, copy vertex locations back
732
static int object_has_edges(Object *ob)
734
if(ob->type==OB_MESH) {
736
if(me->medge) return 1;
738
else if(ob->type==OB_LATTICE) {
746
static void set_body_point(Object *ob, BodyPoint *bp, float *vec)
749
VECCOPY(bp->pos, vec);
750
Mat4MulVecfl(ob->obmat, bp->pos); // yep, sofbody is global coords
751
VECCOPY(bp->origS, bp->pos);
752
VECCOPY(bp->origE, bp->pos);
753
VECCOPY(bp->origT, bp->pos);
755
bp->vec[0]= bp->vec[1]= bp->vec[2]= 0.0;
757
if(ob->softflag & OB_SB_GOAL) {
758
bp->goal= ob->soft->defgoal;
762
/* so this will definily be below SOFTGOALSNAP */
767
bp->contactfrict = 0.0f;
771
/* copy original (new) situation in softbody, as result of matrices or deform */
772
/* is assumed to enter function with ob->soft, but can be without points */
773
static void mesh_update_softbody(Object *ob)
776
MVert *mvert= me->mvert;
777
/* MEdge *medge= me->medge; */ /*unused*/
781
/* possible after a file read... */
782
if(ob->soft->bpoint==NULL) sbObjectToSoftbody(ob);
786
bp= ob->soft->bpoint;
787
for(a=0; a<me->totvert; a++, mvert++, bp++) {
788
VECCOPY(bp->origS, bp->origE);
789
VECCOPY(bp->origE, mvert->co);
790
Mat4MulVecfl(ob->obmat, bp->origE);
791
VECCOPY(bp->origT, bp->origE);
794
if(ob->softflag & OB_SB_EDGES) {
796
/* happens when in UI edges was set */
797
if(ob->soft->bspring==NULL)
798
if(object_has_edges(ob)) sbObjectToSoftbody(ob);
800
/* hrms .. do springs alter their lenght ? (yes, mesh keys would (ton))
802
bs= ob->soft->bspring;
803
bp= ob->soft->bpoint;
804
for(a=0; (a<me->totedge && a < ob->soft->totspring ); a++, medge++, bs++) {
805
bs->len= VecLenf( (bp+bs->v1)->origE, (bp+bs->v2)->origE);
814
static void get_scalar_from_vertexgroup(Object *ob, int vertID, short groupindex, float *target)
815
/* result 0 on success, else indicates error number
816
-- kind of *inverse* result defintion,
817
-- but this way we can signal error condition to caller
818
-- and yes this function must not be here but in a *vertex group module*
824
/* spot the vert in deform vert list at mesh */
825
if(ob->type==OB_MESH) {
826
if (((Mesh *)ob->data)->dvert) {
827
dv = ((Mesh*)ob->data)->dvert + vertID;
828
/* Lets see if this vert is in the weight group */
829
for (i=0; i<dv->totweight; i++){
830
if (dv->dw[i].def_nr == groupindex){
831
*target= dv->dw[i].weight; /* got it ! */
839
/* makes totally fresh start situation */
840
static void mesh_to_softbody(Object *ob)
844
MVert *mvert= me->mvert;
845
MEdge *medge= me->medge;
851
if (ob->softflag & OB_SB_EDGES) totedge= me->totedge;
854
/* renew ends with ob->soft with points and edges, also checks & makes ob->soft */
855
renew_softbody(ob, me->totvert, totedge);
857
/* we always make body points */
860
goalfac= ABS(sb->maxgoal - sb->mingoal);
862
for(a=me->totvert; a>0; a--, mvert++, bp++) {
864
set_body_point(ob, bp, mvert->co);
866
/* get scalar values needed *per vertex* from vertex group functions,
867
so we can *paint* them nicly ..
868
they are normalized [0.0..1.0] so may be we need amplitude for scale
869
which can be done by caller but still .. i'd like it to go this way
872
if((ob->softflag & OB_SB_GOAL) && sb->vertgroup) {
873
get_scalar_from_vertexgroup(ob, me->totvert - a, sb->vertgroup-1, &bp->goal);
874
// do this always, regardless successfull read from vertex group
875
bp->goal= sb->mingoal + bp->goal*goalfac;
877
/* a little ad hoc changing the goal control to be less *sharp* */
878
bp->goal = (float)pow(bp->goal, 4.0f);
880
/* to proove the concept
881
this would enable per vertex *mass painting*
882
strcpy(name,"SOFTMASS");
883
error = get_scalar_from_named_vertexgroup(ob,name,me->totvert - a,&temp);
884
if (!error) bp->mass = temp * ob->rangeofmass;
888
/* but we only optionally add body edge springs */
889
if (ob->softflag & OB_SB_EDGES) {
893
for(a=me->totedge; a>0; a--, medge++, bs++) {
897
bs->len= VecLenf( (bp+bs->v1)->origS, (bp+bs->v2)->origS);
901
/* insert *diagonal* springs in quads if desired */
902
if (ob->softflag & OB_SB_QUADS) {
903
add_mesh_quad_diag_springs(ob);
906
build_bps_springlist(ob); /* big mesh optimization */
912
/* copies current sofbody position in mesh, so do this within modifier stacks! */
913
static void softbody_to_mesh(Object *ob)
920
bp= ob->soft->bpoint;
922
for(a=me->totvert; a>0; a--, mvert++, bp++) {
923
VECCOPY(mvert->co, bp->pos);
924
Mat4MulVecfl(ob->imat, mvert->co); // softbody is in global coords
929
/* makes totally fresh start situation */
930
static void lattice_to_softbody(Object *ob)
933
Lattice *lt= ob->data;
938
totvert= lt->pntsu*lt->pntsv*lt->pntsw;
940
/* renew ends with ob->soft with points and edges, also checks & makes ob->soft */
941
renew_softbody(ob, totvert, 0);
943
/* we always make body points */
946
for(a= totvert, bp= lt->def, bop= sb->bpoint; a>0; a--, bp++, bop++) {
947
set_body_point(ob, bop, bp->vec);
951
/* copies current sofbody position */
952
static void softbody_to_lattice(Object *ob)
955
Lattice *lt= ob->data;
960
totvert= lt->pntsu*lt->pntsv*lt->pntsw;
963
for(a= totvert, bp= lt->def, bop= sb->bpoint; a>0; a--, bp++, bop++) {
964
VECCOPY(bp->vec, bop->pos);
965
Mat4MulVecfl(ob->imat, bp->vec); // softbody is in global coords
969
/* copy original (new) situation in softbody, as result of matrices or deform */
970
/* is assumed to enter function with ob->soft, but can be without points */
971
static void lattice_update_softbody(Object *ob)
973
Lattice *lt= ob->data;
978
totvert= lt->pntsu*lt->pntsv*lt->pntsw;
980
/* possible after a file read... */
981
if(ob->soft->bpoint==NULL) sbObjectToSoftbody(ob);
983
for(a= totvert, bp= lt->def, bop= ob->soft->bpoint; a>0; a--, bp++, bop++) {
984
VECCOPY(bop->origS, bop->origE);
985
VECCOPY(bop->origE, bp->vec);
986
Mat4MulVecfl(ob->obmat, bop->origE);
987
VECCOPY(bop->origT, bop->origE);
992
/* copies softbody result back in object */
993
/* only used in sbObjectStep() */
994
static void softbody_to_object(Object *ob)
997
if(ob->soft==NULL) return;
999
/* inverse matrix is not uptodate... */
1000
Mat4Invert(ob->imat, ob->obmat);
1004
softbody_to_mesh(ob);
1007
softbody_to_lattice(ob);
1012
/* copy original (new) situation in softbody, as result of matrices or deform */
1013
/* used in sbObjectStep() and sbObjectReset() */
1014
/* assumes to have ob->soft, but can be entered without points */
1015
static void object_update_softbody(Object *ob)
1020
mesh_update_softbody(ob);
1023
lattice_update_softbody(ob);
1029
/* return 1 if succesfully baked and applied step */
1030
static int softbody_baked_step(Object *ob, float framenr)
1032
SoftBody *sb= ob->soft;
1033
SBVertex *key0, *key1, *key2, *key3;
1035
float data[4], sfra, efra, cfra, dfra, fac; // start, end, current, delta
1038
/* precondition check */
1039
if(sb==NULL || sb->keys==NULL || sb->totkey==0) return 0;
1040
/* so we got keys, but no bodypoints... even without simul we need it for the bake */
1041
if(sb->bpoint==NULL) sb->bpoint= MEM_callocN( sb->totpoint*sizeof(BodyPoint), "bodypoint");
1043
/* convert cfra time to system time */
1044
sfra= (float)sb->sfra;
1045
cfra= bsystem_time(ob, NULL, framenr, 0.0);
1046
efra= (float)sb->efra;
1047
dfra= (float)sb->interval;
1049
/* offset in keys array */
1050
ofs1= floor( (cfra-sfra)/dfra );
1053
key0=key1=key2=key3= *sb->keys;
1055
else if(ofs1 >= sb->totkey-1) {
1056
key0=key1=key2=key3= *(sb->keys+sb->totkey-1);
1059
key1= *(sb->keys+ofs1);
1060
key2= *(sb->keys+ofs1+1);
1062
if(ofs1>0) key0= *(sb->keys+ofs1-1);
1065
if(ofs1<sb->totkey-2) key3= *(sb->keys+ofs1+2);
1069
sb->ctime= cfra; // needed?
1072
fac= ((cfra-sfra)/dfra) - (float)ofs1;
1073
CLAMP(fac, 0.0, 1.0);
1074
set_four_ipo(fac, data, KEY_BSPLINE);
1076
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++, key0++, key1++, key2++, key3++) {
1077
bp->pos[0]= data[0]*key0->vec[0] + data[1]*key1->vec[0] + data[2]*key2->vec[0] + data[3]*key3->vec[0];
1078
bp->pos[1]= data[0]*key0->vec[1] + data[1]*key1->vec[1] + data[2]*key2->vec[1] + data[3]*key3->vec[1];
1079
bp->pos[2]= data[0]*key0->vec[2] + data[1]*key1->vec[2] + data[2]*key2->vec[2] + data[3]*key3->vec[2];
1082
softbody_to_object(ob);
1087
/* only gets called after succesfully doing softbody_step */
1088
/* already checked for OB_SB_BAKE flag */
1089
static void softbody_baked_add(Object *ob, float framenr)
1091
SoftBody *sb= ob->soft;
1094
float sfra, efra, cfra, dfra, fac1; // start, end, current, delta
1097
/* convert cfra time to system time */
1098
sfra= (float)sb->sfra;
1099
cfra= bsystem_time(ob, NULL, framenr, 0.0);
1100
efra= (float)sb->efra;
1101
dfra= (float)sb->interval;
1104
if(sb->sfra >= sb->efra) return; // safety, UI or py setting allows
1105
if(sb->interval<1) sb->interval= 1; // just be sure
1107
sb->totkey= 1 + (int)(ceil( (efra-sfra)/dfra ) );
1108
sb->keys= MEM_callocN( sizeof(void *)*sb->totkey, "sb keys");
1111
/* now find out if we have to store a key */
1113
/* offset in keys array */
1119
ofs1= floor( (cfra-sfra)/dfra );
1120
fac1= ((cfra-sfra)/dfra) - (float)ofs1;
1122
if( fac1 < 1.0/dfra ) {
1124
key= *(sb->keys+ofs1);
1126
*(sb->keys+ofs1)= key= MEM_mallocN(sb->totpoint*sizeof(SBVertex), "softbody key");
1128
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++, key++) {
1129
VECCOPY(key->vec, bp->pos);
1135
/* ************ Object level, exported functions *************** */
1137
/* allocates and initializes general main data */
1138
SoftBody *sbNew(void)
1142
sb= MEM_callocN(sizeof(SoftBody), "softbody");
1144
sb->mediafrict= 0.5;
1147
sb->physics_speed= 1.0;
1150
sb->goalspring= 0.5;
1160
sb->sfra= G.scene->r.sfra;
1161
sb->efra= G.scene->r.efra;
1167
void sbFree(SoftBody *sb)
1169
free_softbody_intern(sb);
1174
/* makes totally fresh start situation */
1175
void sbObjectToSoftbody(Object *ob)
1180
mesh_to_softbody(ob);
1183
lattice_to_softbody(ob);
1187
if(ob->soft) ob->soft->ctime= bsystem_time(ob, NULL, (float)G.scene->r.cfra, 0.0);
1188
ob->softflag &= ~OB_SB_REDO;
1191
/* reset all motion */
1192
void sbObjectReset(Object *ob)
1194
SoftBody *sb= ob->soft;
1198
if(sb==NULL) return;
1199
if(sb->keys && sb->totkey) return; // only as cpu time saver
1201
sb->ctime= bsystem_time(ob, NULL, (float)G.scene->r.cfra, 0.0);
1203
object_update_softbody(ob);
1205
for(a=sb->totpoint, bp= sb->bpoint; a>0; a--, bp++) {
1206
// origS is previous timestep
1207
VECCOPY(bp->origS, bp->origE);
1208
VECCOPY(bp->pos, bp->origE);
1209
VECCOPY(bp->origT, bp->origE);
1210
bp->vec[0]= bp->vec[1]= bp->vec[2]= 0.0f;
1212
// no idea about the Heun stuff! (ton)
1213
VECCOPY(bp->prevpos, bp->pos);
1214
VECCOPY(bp->prevvec, bp->vec);
1215
VECCOPY(bp->prevdx, bp->vec);
1216
VECCOPY(bp->prevdv, bp->vec);
1221
/* simulates one step. framenr is in frames */
1222
/* copies result back to object, displist */
1223
void sbObjectStep(Object *ob, float framenr)
1229
float ctime, forcetime;
1232
/* this variable is set while transform(). with lattices also having a softbody,
1233
it calls lattice_modifier() all the time... has no displist yet. Is temporal
1234
hack which should be resolved with proper depgraph usage + storage of deformed
1235
vertices in lattice (ton) */
1236
if(G.moving) return;
1238
/* baking works with global time */
1239
if(!(ob->softflag & OB_SB_BAKEDO) )
1240
if(softbody_baked_step(ob, framenr) ) return;
1242
/* remake softbody if: */
1243
if( (ob->softflag & OB_SB_REDO) || // signal after weightpainting
1244
(ob->soft==NULL) || // just to be nice we allow full init
1245
(ob->soft->bpoint==NULL) ) // after reading new file, or acceptable as signal to refresh
1246
sbObjectToSoftbody(ob);
1250
/* still no points? go away */
1251
if(sb->totpoint==0) return;
1253
/* reset deflector cache, sumohandle is free, but its still sorta abuse... (ton) */
1254
for(base= G.scene->base.first; base; base= base->next) {
1255
base->object->sumohandle= NULL;
1258
/* checking time: */
1259
ctime= bsystem_time(ob, NULL, framenr, 0.0);
1260
dtime= ctime - sb->ctime;
1261
// bail out for negative or for large steps
1262
if(dtime<0.0 || dtime >= 9.9*G.scene->r.framelen) { // G.scene->r.framelen corrects for frame-mapping, so this is actually 10 frames for UI
1269
if(dtime > 0.0) { // note: what does this mean now? (ton)
1271
//dtime is still in [frames]
1272
//we made sure dtime is >= 0.0
1273
//but still need to handle dtime == 0.0 -> just return sb as is, just to be nice
1274
object_update_softbody(ob);
1276
if (TRUE) { // RSOL1 always true now (ton)
1277
/* special case of 2nd order Runge-Kutta type AKA Heun */
1278
float timedone =0.0; // how far did we get without violating error condition
1279
/* loops = counter for emergency brake
1280
* we don't want to lock up the system if physics fail
1283
SoftHeunTol = sb->rklimit; // humm .. this should be calculated from sb parameters and sizes
1285
forcetime = dtime; /* hope for integrating in one step */
1286
while ( (ABS(timedone) < ABS(dtime)) && (loops < 2000) )
1288
if (ABS(dtime) > 9.0 ){
1289
if(G.f & G_DEBUG) printf("SB_STEPSIZE \n");
1290
break; // sorry but i must assume goal movement can't be interpolated any more
1293
interpolate_exciter(ob,200,(int)(200.0*(timedone/dtime)));
1294
// do predictive euler step
1295
softbody_calc_forces(ob, forcetime);
1296
softbody_apply_forces(ob, forcetime, 1, NULL);
1297
// crop new slope values to do averaged slope step
1298
softbody_calc_forces(ob, forcetime);
1299
softbody_apply_forces(ob, forcetime, 2, &err);
1300
softbody_apply_goalsnap(ob);
1302
if (err > SoftHeunTol){ // error needs to be scaled to some quantity
1303
softbody_restore_prev_step(ob);
1307
float newtime = forcetime * 1.1f; // hope for 1.1 times better conditions in next step
1309
if (err > SoftHeunTol/2.0){ // stay with this stepsize unless err really small
1310
newtime = forcetime;
1312
timedone += forcetime;
1313
if (forcetime > 0.0)
1314
forcetime = MIN2(dtime - timedone,newtime);
1316
forcetime = MAX2(dtime - timedone,newtime);
1320
// move snapped to final position
1321
interpolate_exciter(ob, 2, 2);
1322
softbody_apply_goalsnap(ob);
1325
if (loops > HEUNWARNLIMIT) /* monitor high loop counts say 1000 after testing */
1326
printf("%d heun integration loops/frame \n",loops);
1330
/* do brute force explicit euler */
1331
/* inner intagration loop */
1333
// loop n times so that n*h = duration of one frame := 1
1334
// x(t+h) = x(t) + h*v(t);
1335
// v(t+h) = v(t) + h*f(x(t),t);
1336
timescale = (int)(sb->rklimit * ABS(dtime));
1337
for(t=1 ; t <= timescale; t++) {
1338
if (ABS(dtime) > 15 ) break;
1340
/* the *goal* mesh must use the n*h timing too !
1341
use *cheap* linear intepolation for that */
1342
interpolate_exciter(ob,timescale,t);
1343
if (timescale > 0 ) {
1344
forcetime = dtime/timescale;
1346
/* does not fit the concept sloving ODEs :) */
1347
/* softbody_apply_goal(ob,forcetime ); */
1349
/* explicit Euler integration */
1350
/* we are not controling a nuclear power plant!
1351
so rought *almost* physical behaviour is acceptable.
1352
in cases of *mild* stiffnes cranking up timscale -> decreasing stepsize *h*
1353
avoids instability */
1354
softbody_calc_forces(ob,forcetime);
1355
softbody_apply_forces(ob,forcetime,0, NULL);
1356
softbody_apply_goalsnap(ob);
1359
/* ok here comes the �berhammer
1360
use a semi implicit euler integration to tackle *all* stiff conditions
1361
but i doubt the cost/benifit holds for most of the cases
1369
/* and apply to vertices */
1370
softbody_to_object(ob);
1373
} // if(ABS(dtime) > 0.0)
1375
// rule : you have asked for the current state of the softobject
1376
// since dtime= ctime - ob->soft->ctime== 0.0;
1377
// and we were not notifified about any other time changes
1379
softbody_to_object(ob);
1382
/* reset deflector cache */
1383
for(base= G.scene->base.first; base; base= base->next) {
1384
if(base->object->sumohandle) {
1385
MEM_freeN(base->object->sumohandle);
1386
base->object->sumohandle= NULL;
1390
if(ob->softflag & OB_SB_BAKEDO) softbody_baked_add(ob, framenr);
added
removed
source/blender/blenkernel/intern/softbody.c
