; Welcome to the FSH mangle sourcecode!
;
; This source code is released under the 'Gee whizz, Batman' license
; If you use this code in a commercial project, you're... you're...
; 'kay, you're just nuts. Check out the Paul Bourke reference I 
; keep throwing around!
;
; If you desperately want to use this stuff, then just give me a smily
; face in the credits and e-mail me about it!
;
; Okie doke, here's the skinny about the world of particles... let's
; see, where to begin?
;
; Particles:
;   * They have a position                 - They've gotta BE somewhere.
;   * They have a mass                     - This is their inertia...
;                                            unless you're going relativistic.
;   * They have a velocity                 - Position / time = velocity. 
;   * They have an acceleration            - Velocity / time = acceleration.
;   
; Springs:
;   * They have a rest length              - Springs want to be this length.
;   * They have a spring constant          - This represents their strength.
;   * This have a damping constant         - This represents their action under
;                                            different environments (I think)
;   * They connect two particles together  - It's important that the springs
;                                            refer to two particles, as opposed
;                                            to the other way round. It's harder
;                                            that way, you can't do Spring force
;                                            on B = - Spring force on A like I
;                                            do here.
;
; Now listen close, 'cause this is the way it's going to work...
;
; Imagine this: You've been given a basketball and a machine that makes time
; completely stop.
; Now, imagine that you hold the basketball upto about waist height, and then
; stop time. the basketball simply hangs there in midair. Now try pushing the thing...
; It doesn't move. Try pushing the ball from another angle. And again. And again.
; Imagine that there's something about the basketball that remembers every direction
; you've pushed and how hard you've pushed, and when you start time again, the ball
; travels out in that direction. 
;
; In this particle system, that's the best way I could come up with to explain
; 'force accumulators'. If you've got a box, and you and a friend push it from
; opposite sides, it won't move. It's like that.
;
; A particle which is connected to a bunch of springs will experience a spring
; action force from each of them. In my system here, every single spring action
; force is added to the particles force accumulator, which is then used to
; calculate the final change in velocity for the particle.
;
; In this system:
;  * The particles and springs are set up in their starting positions.
;  Every frame:
;    * Every particles force accumulator is set to zero.
;    * Every spring has it's deviation from it's rest length
;       computed, and that is used to determine the force it exerts.
;       If it acts on a non-fixed particle, it will add it's force
;       onto the accumulator of the particles it is attached to.
;    * Here, the user can influence the particles by adding their own
;       force onto the accumulators.
;    * Using the final values in the accumulator, acceleration for this frame
;       is calculated using acceleration = force / mass.
;    * Then, velocity change this frame is calculated as 
;       change in velocity this frame = acceleration * time.
;       This is then added onto the particles current velocity.
;    * Then, the change in position per frame is calculated using
;       change in position this frame = velocity * time.
;    * This is then added onto the particles velocity.
;    * Now that the particles have had their motion this frame calculated
;       and applied, you can now render the frame however you like. Remember
;       the internal particle coords don't have to correspond to the BlitzBasic
;       screen coords. (In this case, I've flipped the Y axis for the particle
;       coords. This means that a negative Y in the particle coords means that it's
;       below the screen. (Which means it has a highly positive BlitzBasic screen
;       y coordinate)
;
; Notice that there is a Z dimension calculation here. If you take out everything to
; do with the Z, then this simulation will work just the same. All distances along the
; z axis are defined as zero for this simulation... so they shouldn't change, so there's
; no change in the Z position due to spring action forces (which are directly
; proportional to the vector along which the spring lies... with a zero Z coordinate,
; these should all be parallel to the plane of the screen, so no 'inward' or 'outward'
; motion.
;
; The spring equation is based off hooke's spring law, the explanation of which
; is on the Paul Bourke site! I based this code off of that equation, referring to the
; C library code after I'd wrote it. As it turned out, my code was pretty similar, so I
; guess I got it right! Wahoo!
;
; How does particles and springs make a floating space handkerchief, you ask me?
; 
; Most modelling situations involve the modelling of extremely complex systems by
; observing the interactions of many much smaller systems. It just happens that a
; collection of masses attempting to retain a fixed distance from one another mimics
; the motion of fabric. I guess if you had a bunch of springs and one of these neat
; magenetic ball sets, you could make a real life FSH to play around with.
;
; Exercise: Make your own handkerchief my altering the algorithms at the bottom of
; this file.
;
; Exercise: Alter the algorithm to make a non-rectilinear handkerchief.
;
; Exercise: Collide particles with the edge of the screen.
;
; Exercise: Make it so that springs which exceed a certain length 'snap' and
; are removed from the simulation, causing their attached particles to
; become freed.
;
; Exercise: Call your new game HSL.
;
; Exercise: Rewrite the algorithm to accept particle specifications from an external file
; so you can write simple frame 'machines' in notepad and import them.
;
; Exercise: Make a stickman out of little rectangular boxes joined by common spring
; points at the joints, and float him about.
;
; Exercise: -turvat. Need I say more?
;
; Hope this helps someone! ^_^ Have fun now!
;
; Mathew Carr - mrdictionary.net (matt@mrdictionary.net)

Const width=320,height=240
Const gravity#=0.0981
Const ge=0,j=20                                   ; Gravity enabled? Distance between particles.
Const uniformmass#=0,uniforms#=0,uniformd#=0      ; Set these to something other than zero to change the
                                                  ; starting properties of the springs and particles.

Graphics width,height,16,2
SetBuffer BackBuffer()

SeedRnd MilliSecs()^2

; Here's the physical properties of the particles
; No types here... didn't know how to use them at this point
; I wrote this pre-hypercalm, y'see. I guess types would help
; but you'd have be real good at using the object ID's with
; the spring A and B properties.
Const particles=100
Dim mass#(particles)
Dim xpos#(particles)
Dim ypos#(particles)
Dim zpos#(particles)
Dim xvel#(particles)
Dim yvel#(particles)
Dim zvel#(particles)
Dim xacc#(particles)
Dim yacc#(particles)
Dim zacc#(particles)
Dim xforce#(particles)
Dim yforce#(particles)
Dim zforce#(particles)
Dim fix(particles)
Dim particleexisting(particles)

; Here's the physical properties of the springs
Const springs=500
Dim A(springs)
Dim B(springs)
Dim sconst#(springs)
Dim dconst#(springs)
Dim rest#(springs)
Dim springexisting(springs)

;lets create the particles, and then the springs.

Gosub mathparticle

Gosub mathspring

For cold=0 To particles
	If uniformmass#<>0 Then mass#(cold)=uniformmass#
Next
For cold=0 To springs
	If uniforms#<>0 Then sconst#(cold)=uniforms#
	If uniformd#<>0 Then dconst#(cold)=uniformd#
Next

While Not KeyHit(1)
	ot=MilliSecs()
	Gosub forces
	If MouseDown(1)=0
	If KeyDown(200)=1 And fix(0)=0 Then yforce#(0)=yforce#(0)+1
	If KeyDown(203)=1 And fix(0)=0 Then xforce#(0)=xforce#(0)-1
	If KeyDown(205)=1 And fix(0)=0 Then xforce#(0)=xforce#(0)+1
	If KeyDown(208)=1 And fix(0)=0 Then yforce#(0)=yforce#(0)-1
	Else
	If KeyDown(200)=1 And fix(99)=0 Then yforce#(99)=yforce#(99)+1
	If KeyDown(203)=1 And fix(99)=0 Then xforce#(99)=xforce#(99)-1
	If KeyDown(205)=1 And fix(99)=0 Then xforce#(99)=xforce#(99)+1
	If KeyDown(208)=1 And fix(99)=0 Then yforce#(99)=yforce#(99)-1
	EndIf

	If KeyHit(45)=1
		For doh=0 To dataparticles-1
			xvel#(doh)=0
			yvel#(doh)=0
		Next
	EndIf
	Gosub particles
	Gosub drawlines
	Gosub drawparticles
	Gosub sync
Wend
End

; Call this subroutine to get the system to blank the particles force accumulators
; and fill them with the result of the spring action from every spring.
.forces
; Reset the accumulators on the particles to zero.
; The accumulators will the total all the calculated forces due to the spring action.

	For s=0 To dataparticles-1
		xforce#(s)=0
		yforce#(s)=0
		zforce#(s)=0
	Next
	
; This code will work out the forces due to the spring action for all springs.
	
	For s=0 To datasprings-1
		If springexisting(s)=1
			start=A(s)
			finish=B(s)
			dx#=(xpos#(start)-xpos#(finish))
			dy#=(ypos#(start)-ypos#(finish))
			dz#=(zpos#(start)-zpos#(finish))
			length#=Sqr(dx#^2+dy#^2+dz#^2)+0.000000000001 ; No div/0
			                                              ; This will happen if two particles
			                                              ; which are connected happen to be
			                                              ; at the same point simultaneously.
			fx#=4*sconst#(s)*(length#-rest#(s))
			fx#=fx#+dconst#(s)*(xvel#(start)-xvel#(finish))*dx#/length#
			fx#=fx#*(0-dx#)/length#
			fy#=4*sconst#(s)*(length#-rest#(s))
			fy#=fy#+dconst#(s)*(yvel#(start)-yvel#(finish))*dy#/length#
			fy#=fy#*(0-dy#)/length#
			fz#=4*sconst#(s)*(length#-rest#(s))
			fz#=fz#+dconst#(s)*(zvel#(start)-zvel#(finish))*dz#/length#
			fz#=fz#*(0-dz#)/length#
			If fix(start)=0                               ; If the start is fixed, don't accumulate.
				xforce#(start)=xforce#(start)+fx#         
				yforce#(start)=yforce#(start)+fy#         ; These three add the calculated forces
				zforce#(start)=zforce#(start)+fz#         ; to the accumulators of the particles on the spring.
			EndIf
			If fix(finish)=0                              ; If the end is fixed, don't accumulate.
				xforce#(finish)=xforce#(finish)-fx#
				yforce#(finish)=yforce#(finish)-fy#       ; These three add the calculated forces again
				zforce#(finish)=zforce#(finish)-fz#       ; Force on particle B is equal and opposite to the
			EndIf                                         ; force on particle A.
		EndIf
	Next
	If ge=1                                               ; If gravity is enabled, and the particle isn't fixed....
		For s=0 To dataparticles-1
			If particleexisting(s)=1 And fix(s)=0 Then yforce#(s)=yforce#(s)-gravity#*mass#(s)
		Next                                              ; Add gravity to the force accumulator of the particle.
	EndIf
Return


; Call this subroutine to get the system to use the force accumulators to move all the particles.
.particles
	;The acceleration of this frame is the force/mass, as well you should know =P
	For p=0 To dataparticles-1
		xacc#(p)=xforce#(p)/mass#(p)
		yacc#(p)=yforce#(p)/mass#(p)
		zacc#(p)=zforce#(p)/mass#(p)
		xvel#(p)=xvel#(p)+xacc#(p)     ; Here, we should have actually multiplied the acc by time, as
		yvel#(p)=yvel#(p)+yacc#(p)     ; acceleration * time = change in velocity. But, since these are
		zvel#(p)=zvel#(p)+zacc#(p)     ; all inherently unitless anyway, we can assume that the acceleration
		                               ; can be pixels per frame squared, and therefore velocity is pixels
		                               ; per frame, and since we're assuming that this is '1 frame' of motion
		                               ; it all works out.		
		
		xpos#(p)=xpos#(p)+xvel#(p)     ; As above, change in position should equal velocity multiplied by time.
		ypos#(p)=ypos#(p)+yvel#(p)     ; But the time in this case is again '1 frame'.
		zpos#(p)=zpos#(p)+zvel#(p)
	Next
Return

; Call this function to use the BlitzBasic line function to draw all the springs as lines
; don't worry about the calculations... it simply changes from particle coords
; to screen coords (you could use matrices here, if you were feeling clever and
; wanted to realtime rotozoom the particles) and then it checks to see if the 
; length is less or greater than the restlength. There are two colouring algorithms
; which colour it red as the springs get compressed, and green as they get
; extended. 
.drawlines
	For p=0 To datasprings-1
		If springexisting(p)=1
			start=A(p)
			finish=B(p)
			xda=xpos#(start)
			xdb=xpos#(finish)
			yda=height-ypos#(start)
			ydb=height-ypos#(finish)
			dx#=(xpos#(start)-xpos#(finish))
			dy#=(ypos#(start)-ypos#(finish))
			dz#=(zpos#(start)-zpos#(finish))
			length#=Sqr(dx#^2+dy#^2+dz#^2)+0.000000000001
			stress#=length#/rest#(p)
			If stress#<1
			o#=tehcrunchnessconst#
			st#=stress#
			rogerx#=(st#-o#)/(1-o#)
			rad#=255
			gad#=rogerx#*510-255
			bad#=0
			EndIf
			If stress#>1
			o#=tehbreaknessconst#
			st#=stress#
			rogerx#=((st#-1)/(o#-1))
			rad#=255-(rogerx#*255*2)
			gad#=255
			bad#=0
			EndIf
			Color rad,gad,bad
			Line xda,yda,xdb,ydb
		EndIf
	Next
Return

.drawparticles
	For p=0 To dataparticles-1
		Color 255,255,255
	    Oval xpos#(p),height-ypos#(p),2,2,True
	Next
Return

.Sync
    Color 255,255,255
	Text 0,0,1000/(MilliSecs()-ot)
	Flip
	Cls
Return

; This subroutine generate a lattice of particles, by setting their physical properties
; such as mass and position.
.mathparticle
	dataparticles=100
	For yplace=0 To 9
		For xplace=0 To 9
			pid=xplace+yplace*10
			xpos#(pid)=j+j*xplace
			ypos#(pid)=height-(j+j*yplace)
			zpos#(pid)=0
			xvel#(pid)=0
			yvel#(pid)=0
			zvel#(pid)=0
			xforce#(pid)=0
			yforce#(pid)=0
			zforce#(pid)=0
			fix(pid)=0
			mass#(pid)=1
			xacc#(pid)=0
			yacc#(pid)=0
			zacc#(pid)=0
			particleexisting(pid)=1
		Next
	Next
Return

; This subroutine generates a lattice connecting the lattice
; of particles described above. It's a little crap... there's
; a much better one in fsh3d... but this'll work for now.
.mathspring
	datasprings=500
	For yplace=0 To 8
		For xplace=0 To 8
			sid=4*(xplace+yplace*10)
			pid=xplace+yplace*10
			A(sid)=pid
			B(sid)=pid+1
			sconst#(sid)=0.01
			dconst#(sid)=0.01
			springexisting(sid)=1
			rest#(sid)=j
			
			sid=sid+1
			A(sid)=pid
			B(sid)=pid+10
			sconst#(sid)=0.01
			dconst#(sid)=0.01
			springexisting(sid)=1
			rest#(sid)=j
			
			sid=sid+1
			A(sid)=pid+10
			B(sid)=pid+1
			sconst#(sid)=0.01
			dconst#(sid)=0.01
			springexisting(sid)=1
			rest#(sid)=Sqr(j^2+j^2)
			
			sid=sid+1
			A(sid)=pid
			B(sid)=pid+11
			sconst#(sid)=0.01
			dconst#(sid)=0.01
			springexisting(sid)=1
			rest#(sid)=Sqr(j^2+j^2)
		Next
	Next
	A(450)=90 : B(450)=91 : sconst#(450)=0.01 : dconst#(450)=0.01
	springexisting(450)=1 : rest#(450)=j
	A(451)=91 : B(451)=92 : sconst#(451)=0.01 : dconst#(451)=0.01
	springexisting(451)=1 : rest#(451)=j
	A(452)=92 : B(452)=93 : sconst#(452)=0.01 : dconst#(452)=0.01
	springexisting(452)=1 : rest#(452)=j
	A(453)=93 : B(453)=94 : sconst#(453)=0.01 : dconst#(453)=0.01
	springexisting(453)=1 : rest#(453)=j
	A(454)=94 : B(454)=95 : sconst#(454)=0.01 : dconst#(454)=0.01
	springexisting(454)=1 : rest#(454)=j
	A(455)=95 : B(455)=96 : sconst#(455)=0.01 : dconst#(455)=0.01
	springexisting(455)=1 : rest#(455)=j
	A(456)=96 : B(456)=97 : sconst#(456)=0.01 : dconst#(456)=0.01
	springexisting(456)=1 : rest#(456)=j
	A(457)=97 : B(457)=98 : sconst#(457)=0.01 : dconst#(457)=0.01
	springexisting(457)=1 : rest#(457)=j
	A(458)=98 : B(458)=99 : sconst#(458)=0.01 : dconst#(458)=0.01
	springexisting(458)=1 : rest#(458)=j
	A(459)=9 : B(459)=19 : sconst#(459)=0.01 : dconst#(459)=0.01
	springexisting(459)=1 : rest#(459)=j
	A(460)=19 : B(460)=29 : sconst#(460)=0.01 : dconst#(460)=0.01
	springexisting(460)=1 : rest#(460)=j
	A(461)=29 : B(461)=39 : sconst#(461)=0.01 : dconst#(461)=0.01
	springexisting(461)=1 : rest#(461)=j
	A(462)=39 : B(462)=49 : sconst#(462)=0.01 : dconst#(462)=0.01
	springexisting(462)=1 : rest#(462)=j
	A(463)=49 : B(463)=59 : sconst#(463)=0.01 : dconst#(463)=0.01
	springexisting(463)=1 : rest#(463)=j
	A(464)=59 : B(464)=69 : sconst#(464)=0.01 : dconst#(464)=0.01
	springexisting(464)=1 : rest#(464)=j
	A(465)=69 : B(465)=79 : sconst#(465)=0.01 : dconst#(465)=0.01
	springexisting(465)=1 : rest#(465)=j
	A(466)=79 : B(466)=89 : sconst#(466)=0.01 : dconst#(466)=0.01
	springexisting(466)=1 : rest#(466)=j
	A(467)=89 : B(467)=99 : sconst#(467)=0.01 : dconst#(467)=0.01
	springexisting(467)=1 : rest#(467)=j
Return
; The reason this lot have to be done manually is that this algorithm    -----
; treats the lattice as a bunch of squares with this sort of shape:     |\   /
;                                                                       | \ /
; The right hand edge and the bottom edge of the last row and column    |  X
; have to specified manually. I fixed this in FSH3d though.             | / \
;                                                                       |/   \

;

End