3D Graphics engine (Aug 28 update): semi-intelligent wireframe plots

[STxAxTIC]

(UPDATED Aug 28 2011) Final summer update: progress on intelligent wireframes.

(UPDATED Aug 02 2011) Added some non-trivial functionality to be useful in scientific graphing and research. Check out option "9" for the "single slit experiment on a rainy day".

(UPDATED July 29 2011) I've added a "water engine" capability (inspired by DarthWho) with nonreflecting boundary conditions, which was actually tricky.

(UPDATED July 25 2011) This might be the final version of my 3D graphics... no longer a calculator... engine.

This copy has the full functionality of plotting full 3d molecules, curves, and surfaces. It utilizes time evolution, a relaxation engine, and wave propagation, all done in accordance with Laplace's equation and the classical wave equation. Every last thing in this program is adjustable, and you can change the display when running.

Give it a few tries. Make sure you press T to begin the animations! (It's also worth reading the first few instructions when they pop up.)

For output samples ---> http://people.umass.edu/wfbarnes/3dcollageweb.png
THE MAIN SOURCE CODE IS AVAILABLE HERE---> http://people.umass.edu/wfbarnes/3DTRWGRAPH_FB64.BAS
If you are feeling lazy... ---> http://people.umass.edu/wfbarnes/3DTRWGRAPH_FB64.exe

Use it without my permission if you want. If you understand the code, then that's as good as asking permission.

For anyone interested in learning how this is basically done, I've prepared a simple demo below: (CLICK THE LINK FOR MAIN CODE)

Code: [Select]

SCREEN 12

globaldelay = .5 * 10 ^ 5 'CHANGE THIS NUMBER FOR YOUR MACHINE IF NEEDED

numvectors = 3000 'CHANGE THIS NUMBER FOR YOUR MACHINE IF NEEDED

'*****************************************************************************
'DEFINE WHAT OUR CUBE LOOKS LIKE HERE
'Later in the code, the points we define here will be called "atoms". (Why not?)
'DATA x, y, z, color, point index, neighbor 1, neighbor 2, etc. (up to 5 neighbors)
'I was a little sloppy about being systematic with neighbor connections...
DATA -1.0,-1.0,-1.0,14,1,5,0,0,0,0
DATA -1.0,-1.0,1.00,04,2,1,4,0,0,0
DATA -1.0,1.00,-1.0,09,3,1,4,0,0,0
DATA -1.0,1.00,1.00,08,4,3,0,0,0,0
DATA 1.00,-1.0,-1.0,15,5,7,0,0,0,0
DATA 1.00,-1.0,1.00,07,6,5,2,0,0,0
DATA 1.00,1.00,-1.0,12,7,8,3,0,0,0
DATA 1.00,1.00,1.00,05,8,4,6,0,0,0
'*****************************************************************************

'screen unit vectors in 3-space (uhat, vhat, nhat are not visual)
'uhat is the x-direction for your screen.
'vhat is the y-direction for your screen.
'nhat is an imaginary line that points right at you, out of the screen.
DIM uhat(1 TO 3), vhat(1 TO 3), nhat(1 TO 3)

'main vectors in 3-space (xcoord, ycoord, zcoord, color)
'This array is where you actually will store the information about your object.
'In this case, it will be a cube.
DIM vec(numvectors, 4)

'main vector screen projections in 2-space
'This is where your "flattened" (projected) object is stored.
'Reminder: the object (cube in this case) is flattened onto the screen.
DIM vecp(numvectors, 1 TO 2)
DIM vecp.old(numvectors, 1 TO 2)
DIM vecdotnhat(numvectors)
DIM vecdotnhat.old(numvectors)

'vectors tangent to screen plane (defines initial camera angle)
uhat(1) = -3: uhat(2) = 5: uhat(3) = 1 / 4
vhat(1) = -1: vhat(2) = -1: vhat(3) = 8

COLOR 15: LOCATE 2, 15: PRINT "Use keys W, S, A, D, & Q, E, C, V to rotate display."

zoom = 38.8
speedconst = 25
NumAtoms = 8
numvectors = NumAtoms

REDIM vec(numvectors, 10)

'generate main vectors (xcoord, ycoord, zcoord, color, neighbor 1, neighbor 2, ...)

'This loop will store the cube as an array.
'It stores the location of each corner, and the color of that corner...
'AND it reads which points connect to which neighbors.
FOR i = 1 TO numvectors
    READ vec(i, 1), vec(i, 2), vec(i, 3), vec(i, 4), vec(i, 5), vec(i, 6), vec(i, 7), vec(i, 8), vec(i, 9), vec(i, 10)
NEXT
 
'begin main loop
DO
    igdel = 0: DO: igdel = igdel + .01: LOOP UNTIL igdel > globaldelay
    GOSUB redraw
    GOSUB keyprocess
LOOP

keyprocess:
SELECT CASE LCASE$(INKEY$)
    CASE "d": GOSUB rotate.uhat.plus
    CASE "a": GOSUB rotate.uhat.minus
    CASE "w": GOSUB rotate.vhat.plus
    CASE "s": GOSUB rotate.vhat.minus
    CASE "q": GOSUB rotate.counterclockwise
    CASE "e": GOSUB rotate.clockwise
    CASE "c": GOSUB rotate.uhat.plus: GOSUB normalize.screen.vectors: GOSUB rotate.counterclockwise
    CASE "v": GOSUB rotate.uhat.minus: GOSUB normalize.screen.vectors: GOSUB rotate.clockwise
    CASE CHR$(27): END
END SELECT
RETURN

redraw:
GOSUB normalize.screen.vectors
GOSUB compute.screen.normal
GOSUB projectvectors
GOSUB PlotEverything
GOSUB storeprojections
RETURN

convert:
'convert graphics from cartesian coordinates to screen coordinates
x0 = x: y0 = y
x = x0 + 320
y = -y0 + 240
RETURN

normalize.screen.vectors:
'normalize the two vectors that define the screen orientation
'To "normalize" means to make the vectors have length 1.
uhatmag = SQR(uhat(1) ^ 2 + uhat(2) ^ 2 + uhat(3) ^ 2)
uhat(1) = uhat(1) / uhatmag: uhat(2) = uhat(2) / uhatmag: uhat(3) = uhat(3) / uhatmag
vhatmag = SQR(vhat(1) ^ 2 + vhat(2) ^ 2 + vhat(3) ^ 2)
vhat(1) = vhat(1) / vhatmag: vhat(2) = vhat(2) / vhatmag: vhat(3) = vhat(3) / vhatmag
uhatdotvhat = uhat(1) * vhat(1) + uhat(2) * vhat(2) + uhat(3) * vhat(3)
IF SQR(uhatdotvhat ^ 2) > .001 THEN CLS: LOCATE 5, 5: PRINT "Screen vectors are not perpendicular.": END
RETURN

compute.screen.normal:
'compute nhat, the normal to the screen vectors, using cross product
'I won't even explain this unless you know some vector calculus.
'Just trust it.
nhat(1) = uhat(2) * vhat(3) - uhat(3) * vhat(2)
nhat(2) = uhat(3) * vhat(1) - uhat(1) * vhat(3)
nhat(3) = uhat(1) * vhat(2) - uhat(2) * vhat(1)
nhatmag = SQR(nhat(1) ^ 2 + nhat(2) ^ 2 + nhat(3) ^ 2)
nhat(1) = nhat(1) / nhatmag: nhat(2) = nhat(2) / nhatmag: nhat(3) = nhat(3) / nhatmag
RETURN

projectvectors:
'project main vectors onto the screen plane & compute other quantities
'Again, project means "flatten".
FOR i = 1 TO numvectors
    vecp(i, 1) = vec(i, 1) * uhat(1) + vec(i, 2) * uhat(2) + vec(i, 3) * uhat(3)
    vecp(i, 2) = vec(i, 1) * vhat(1) + vec(i, 2) * vhat(2) + vec(i, 3) * vhat(3)
    vecdotnhat(i) = vec(i, 1) * nhat(1) + vec(i, 2) * nhat(2) + vec(i, 3) * nhat(3)
    'normalize
    magvec = (SQR(vec(i, 1) ^ 2 + vec(i, 2) ^ 2 + vec(i, 3) ^ 2))
    IF magvec = 0 THEN magvec = .0000001
    vecdotnhat(i) = vecdotnhat(i) / magvec
NEXT
RETURN

rotate.uhat.plus:
uhat(1) = nhat(1) + speedconst * uhat(1)
uhat(2) = nhat(2) + speedconst * uhat(2)
uhat(3) = nhat(3) + speedconst * uhat(3)
RETURN

rotate.uhat.minus:
uhat(1) = -nhat(1) + speedconst * uhat(1)
uhat(2) = -nhat(2) + speedconst * uhat(2)
uhat(3) = -nhat(3) + speedconst * uhat(3)
RETURN

rotate.vhat.plus:
vhat(1) = nhat(1) + speedconst * vhat(1)
vhat(2) = nhat(2) + speedconst * vhat(2)
vhat(3) = nhat(3) + speedconst * vhat(3)
RETURN

rotate.vhat.minus:
vhat(1) = -nhat(1) + speedconst * vhat(1)
vhat(2) = -nhat(2) + speedconst * vhat(2)
vhat(3) = -nhat(3) + speedconst * vhat(3)
RETURN

rotate.counterclockwise:
v1 = vhat(1): v2 = vhat(2): v3 = vhat(3)
vhat(1) = uhat(1) + speedconst * vhat(1)
vhat(2) = uhat(2) + speedconst * vhat(2)
vhat(3) = uhat(3) + speedconst * vhat(3)
uhat(1) = -v1 + speedconst * uhat(1)
uhat(2) = -v2 + speedconst * uhat(2)
uhat(3) = -v3 + speedconst * uhat(3)
RETURN

rotate.clockwise:
v1 = vhat(1): v2 = vhat(2): v3 = vhat(3)
vhat(1) = -uhat(1) + speedconst * vhat(1)
vhat(2) = -uhat(2) + speedconst * vhat(2)
vhat(3) = -uhat(3) + speedconst * vhat(3)
uhat(1) = v1 + speedconst * uhat(1)
uhat(2) = v2 + speedconst * uhat(2)
uhat(3) = v3 + speedconst * uhat(3)
RETURN

storeprojections:
xhatp.old(1) = xhatp(1): xhatp.old(2) = xhatp(2)
yhatp.old(1) = yhatp(1): yhatp.old(2) = yhatp(2)
zhatp.old(1) = zhatp(1): zhatp.old(2) = zhatp(2)
FOR i = 1 TO numvectors
    vecp.old(i, 1) = vecp(i, 1)
    vecp.old(i, 2) = vecp(i, 2)
NEXT
RETURN

PlotEverything:
'FOR i = numvectors TO 1 STEP -1
FOR i = 1 TO numvectors
    'erase old lines between atoms
    FOR j = 6 TO 10
        IF vec(i, j) = 0 THEN
            EXIT FOR
        ELSE
            x = zoom * vecp.old(i, 1): y = zoom * vecp.old(i, 2)
            GOSUB convert
            x1 = x: y1 = y
            x = zoom * vecp.old((vec(i, j)), 1): y = zoom * vecp.old((vec(i, j)), 2)
            GOSUB convert
            x2 = x: y2 = y
            LINE (x1, y1)-(x2, y2), 0
        END IF
    NEXT
    'erase old atoms
    x = zoom * vecp.old(i, 1): y = zoom * vecp.old(i, 2)
    GOSUB convert
    IF x > 0 AND x < 640 AND y > 0 AND y < 480 THEN
        CIRCLE (x, y), 3, 0
    END IF
    'draw new lines between atoms
    FOR j = 6 TO 10
        IF vec(i, j) = 0 THEN
            EXIT FOR
        ELSE
            x = zoom * vecp(i, 1): y = zoom * vecp(i, 2)
            GOSUB convert
            x1 = x: y1 = y
            x = zoom * vecp((vec(i, j)), 1): y = zoom * vecp((vec(i, j)), 2)
            GOSUB convert
            x2 = x: y2 = y
            LINE (x1, y1)-(x2, y2), vec(i, 4)
        END IF
    NEXT
    'draw new atoms
    x = zoom * vecp(i, 1): y = zoom * vecp(i, 2)
    GOSUB convert
    IF x > 0 AND x < 640 AND y > 0 AND y < 480 THEN
        CIRCLE (x, y), 3, vec(i, 4)
    END IF
    vecdotnhat.old(i) = vecdotnhat(i)
NEXT
RETURN



[Mrwhy]

That is a clever prog, Stx
Any chance we could increase the wind speed and see a bit of turbulence and chaotic flutter?
Perhaps having the flag hang under gravity would help, as the shape of the flag would then feedback on the wind direction and hence flutter-forces.

[Helmet]

Wow! 

3D with QB64 always blows my mind.

[STxAxTIC]

Hey guys, thanks for replying to this post. I'll reposnd to both questions/comments here.

1) Can we add chaotic flutter and gravity? Absolutely yes, and I'm working on it now. In physicist language, I'm adding support for the mesh overlays to obey Laplaces equation and the 2D Wave equation by means of a relaxation technique, so the movements will be realistic and natural. After this point, adding gravity or other forces and boundary conditions should be basically trivial. Adding a "flagpole" on the left, so that the left edge of the flag remains rigid, is easy as well... We'll see if I have time for doing all this...

2) This is not a QB64 program actually - its pure QBasic code. (Wrong forum, I know... shoot me if I must be shot.)

I posted the code in its entirety so far, so uncommenting the "a$ = ..." line should give you full access to the other graphing modes.

Thanks again fellas! I'll update the post if I find the real time needed to add everything I mentioned. I also need to wory about my "real" work as well, haha.

[Mrwhy]

Stx, I admire your courage. Turbulent flow (of air or indeed anything) has been a huge problem for hundreds of years - handled only by experiment. The classic theory, accepted because nobody knew how to test it, was shown wrong by Roulle and Takins - the instabilities (after the first) do not occur one by one as speed increases.

[DarthWho]

Damn where did I put my P90...  

interesting;

i'll post a link to one of the sites I found that had a tutorial on cloth simulation if I can find the site.

[Mrwhy]

Darth, Cloth is very interesting, for apart from being porous and offering a surface to attach the air boundary-layer, it requires ENERGY to bend it and converts some of it into heat.
The bending deflects the wind - a non-linear feedback mechanism - so  in a STEADY wind the flag might well SELF_ORGANISE itself and the local air-flow! Local stability AMID turbulence! - you can see this on the surface of any river.

[STxAxTIC]

Hey guys,

I figured I'd supply an update on the latest progress. As of my last post, I've made some drastic yet subtle-looking improvements, including support for the Laplace equation, so this engine now has a "landscape generator" quality about it. The next edit will include the full-blown wave equation.

If you run this code, choose option "6" for the newest feature, and press "t" to make it animate. As usual, you can rotate and play with depth controls.

To get back to the American flag, use option "4".

Option "5" makes a pretzel when you press "t", it's fun!

And by all means, explore other options too. Change "globaldelay" to suit your machine.

Next update will be relatively soon. Thanks for discussing this!

Code: [Select]

see below


[Mrwhy]

Graet Stuff, STx.
The wave equation will give you damped nonlinear travelling-waves partly reflected at the flag edges. Damped flag-resonances DRIVEN non-linearly by wind vortices. Mutual control-and-excitation by feedback DISTRIBUTED network-like over the surface of the flag.
Worth a few Nobel prizes!

[STxAxTIC]

SEE TOP POST FOR CURRENT VERSION AND COMMENTS

[Clippy]

Please REFRAIN from putting code in the WRONG FORUM! This is the SECOND TIME!

If EVERYBODY did that, then WHAT would ANY forum title names be good for? BTW, you just got here!

[STxAxTIC]

Hey CLIPPY,

Judge my work before you judge my decision to post here. The other forum gets way less traffic, so I ask, why apply the boundaries you suggest, even when you type in CAPS? You know, italics are the standard for getting a point across.

Don't troll your own forums, either. It's unprofessional.

I think you're just jealous.

[Clippy]

Listen BUCCO, I don't need to try to find your crap in the wrong forum. You knew better but did it anyway!

Put your crap in the right forum!

[STxAxTIC]

Stop being a baby Clippy. If I bother you so much, delete my post, or ask Pete to do it. (sup Pete)

If you can even understand one thing I'm doing in that program besides CLS and PRINT, then maybe you'd appreciate solid work when you see it. I see the junk you make. Even if you post it in the perfect forum, it doesn't make up for the fact that you can't even touch what I do.

Tell ya what, if I write _FULLSCREEN at the top, and make it offically QB64, will that make you happy? Let me know, cause I don't wanna hear you whine anymore. Beyond that, I'm not responding to your pettiness, even if you feel that you own this place. Is ther a forum for trolling? Cause your CRAP belongs there.

[TerryRitchie]

This program is amazing!