Jump to content


Member Since 28 Feb 2013
Offline Last Active May 17 2017 03:40 PM

#5248776 Programming scientific GUI's, data and gui layout?

Posted by cadjunkie on 25 August 2015 - 09:10 AM

10 million isn't that much, 10 million floats is around 80MB, you can store that 10 times in just 1GB.


It seems you might not be familiar with finite element models. It's more than 10 million floats. I'm sure it's double precision, and each node has 6 degrees of freedom. Plus, each node can be tied to multiple elements, so there's element data there along with different types of stress and strain data. It's not uncommon for the output files to be 150+ GB, depending on what information is output.

#5248580 Programming scientific GUI's, data and gui layout?

Posted by cadjunkie on 24 August 2015 - 12:04 PM

You should avoid to store original and derived data into the same object. Treat it like variables in a programming language: You have a variable with the original data, you apply an operator, and yield in a result that is stored in another variable.


Not questioning the soundness of this advice (because I do think it's sound), but would that really be feasible on large data sets? Postprocessing something like finite element result data with 1 million nodes is very standard and larger models with 10+ million nodes are common too. I can't imagine trying to have 2 copies of that data in memory. I would think the original data is stored to disk and only 1 copy is in memory and gets operated on. If need be, then it gets reloaded. But maybe I'm wrong though...it's happened once or twice tongue.png

#5247445 How would I solve this polynomial equation system?

Posted by cadjunkie on 18 August 2015 - 12:50 PM

So there are i,n,m triples for which it cannot be solved, but is there a determined way to find out x,y,z solutions for provided i,n,m, or provide empty solution set ?



The only thing I can think of is to solve for z in the 3rd equation in terms of x using the quadratic equation: \( z = \frac{-mx \pm \sqrt{(m^2-4)x^2+8}}{2} \) and then substitute that into the 2nd equation. Then, subtract that result from the 1st equation*, then solve for x. Then, you can solve for z, and then y. 


*EDIT: Even this way I couldn't eliminate all the y's from the equation. Maybe there's some trick there, but I can't see it.


If you're asking for a general analytical (not numerical) technique, I'm not sure what it would be.

#5235221 Engineering vs Programming?

Posted by cadjunkie on 16 June 2015 - 06:33 PM

I work in aeronautics as a structural engineer building large military UAVs. I can't really speak about software jobs personally but IMHO, aeronautical engineering isn't really as glorious as some people make it out to be. Sure, you get to work on some high-profile projects and use some fun software, but ask anyone in aerospace what they *actually* do and you'll usually hear something like "You know the X aircraft, right? Well on the wing there's a series of smaller Ys that do Z and I analyzed several of those" or "I designed the tooling fixture that locate the Ys on aircraft X". Sometimes you do get to do the "cool stuff", but that's usually not the norm and it takes a long while before anyone will trust you enough (i.e. promote you high enough) to design the important (read: cool) things. I personally like the challenges I get faced with and how much I'm learning, but that's not to say that software is any less interesting or challenging. Plus, unlike what your parents seem to think, aerospace engineers look at computers all day. CAD, simulation software, and Microsoft Office-type work is all we ever do, really (other than the other menial work of looking up specs or reading tech manuals online). Also, engineers aren't really people-persons. I'm an extrovert, but I have to go down to the shop floor to get any real human interaction (okay it's not all that bad, but you get my point). Also, more often than not, you're not going to be picking up a wrench or going into the machine shop and building anything. It's numbers work. That's why you get a higher salary than the shop guys. Plus, IMO, programming can be more stable. Aerospace is about having lots of money to do new things. It takes millions to build a simple airframe, so when the economy's good there's a lot of work to do. When it tanks (like in the US in 2008), there's a lot of unemployment. It's a roller coaster. The upshot is that stuff has to get made even in a down economy, so usually there's something in consumer products you can find work in. I'm not trying to make this sound unattractive, but just paint things a bit more realistically so you can get a better picture of what the choice looks like.


I love to develop software (mostly mechanics simulation software) because of the problems I get to try to solve, but I couldn't ever see myself doing that professionally. IMO, if you make your hobby your job, it can get old fast. Sure, some days you'll love to go to work and get excited about it, but on the days you don't want to do it, you have to force yourself to do it. I like to unwind by thinking about something other than work, which programming lets me do. However, it's obvious that enough people love it more than I do that they've made their career doing it. It's really all about what you'd rather do all day.

#5233935 Aggregation of scaling

Posted by cadjunkie on 09 June 2015 - 04:51 PM

Yes, it's doable and surprisingly easy. Any symmetric matrix can apply non-uniform scaling, where the scale factors are the eigenvalues of the matrix and the associated eigenvectors of the matrix are the vectors along which the scaling is applied. You can keep track of the eigenvectors (the scaling vectors) and make the length of each vector equal to the eigenvalue (the scale factor). Then you can construct the symmetric matrix by constructing the regular diagonal scaling matrix (formed from the eigenvalues along the diagonal) and left-multiplying it by the change of basis matrix (formed from the normalized eigenvectors) and right-multiplying by the COB matrix inverse. The result will be the non-uniform scaling matrix you're looking for. 


These posts might help clarify the math I tried to explain:





Hope that helps!

#5225112 Better way for geometry problem

Posted by cadjunkie on 23 April 2015 - 12:44 PM

The easiest way I can think of is to use vector math to find BC (BA+AC=BC) and AD (AD = DB-BA) and get the intersection E from them. You can find the analogue to point C on the circle around D (call it G) to construct a line BG and F would be the intersection of BG and AD. 


I wrote an article that details a method that simply uses the cross product to determine the intersection of 2 lines. The lines need to be in an implicit form (ax+by+c=0) and the resulting point is in homogeneous coordinates (x*w,y*w,w), but those transformations are very simple. This method also doesn't yield any infinities if the lines are parallel, and no trigonometric functions are used.



#5223446 Can i do this using the fourier transform?

Posted by cadjunkie on 15 April 2015 - 09:57 AM

I'm unsure on what you're trying to accomplish. Can you provide a more detailed and concrete example?

#5222938 Understanding Generalized Barycentric Coordinates

Posted by cadjunkie on 13 April 2015 - 09:09 AM

Randy Gaul is right. You lose the one-to-one mapping because you can find other sets of weights that satisfy all the equations.


However, if I understand the application right (which maybe I don't), I don't think you need to care about that for a rasterizer. Losing the one-to-one mapping in this case means that you could have two points be rendered as the same color, but does that matter? I wouldn't think so, but then I might be missing the point of the application. 

#5222136 Understanding Generalized Barycentric Coordinates

Posted by cadjunkie on 08 April 2015 - 04:54 PM

Your example doesn't work because of 2 things:

  1. You're using the normal Cartesian coordinates of the triangle points instead of converting them to barycentric coordinates.
  2. You're basically pre-selecting the barycentric coordinates. You need to use the formulas to compute all the \( w_j \) values and then you'll see \( \sum_{j=0}^n w_j(\mathbf{q}_j-\mathbf{p}) = 0 \), where \( \mathbf{q}_j, \mathbf{p} \) are in Cartesian coordinates.  

In order to get the proper interpolation you're looking for, You have to use the formulas in that paper to convert the point to the barycentric coordinates (the \(\alpha_j\) vector) and then dot that with the values you want those points to have in order to get the interpolated value of the given point.


I will agree with Randy Gaul that you do lose coordinate uniqueness if the polytope isn't a simplex. However, I disagree with him that this paper is not interesting. True, barycentric coordinates with a simplex will guarantee a unique coordinate value for any point, but there are a lot of uses beyond just finding points in a simplex, although you have to be sufficiently interested in computational geometry to appreciate them. In a game programming setting, there might be only a couple of relevant uses beyond this one. 

#5219365 Cubic interpolation on an irregular grid (triangle mesh)?

Posted by cadjunkie on 26 March 2015 - 11:06 AM

As a quick note about the Bezier triangle, you have to be somewhat careful in "ghosting" your vertices. Interpolation is a way of getting continuous data out of a smaller set of discrete data. You've got a certain amount of information in your data set and interpolation won't increase that amount of information. The linear, quadratic, and cubic interpolation schemes are just ways of guessing what's happening between the data points, but it doesn't give you any new information that you didn't have before. I guess the main point here is if you're trying to linearly interpolate the triangle vertices to get these "ghost" points and then use a cubic interpolation to get the colors, you're going to end up with a cubic approximation of a linear interpolation.


The polyharmonic spline is probably a good way to go if you're comfortable with the approximations it's making. It assumes that the grid points are exact and that it minimizes the "bending energy". The way Eberly puts it, it's like a thin metal sheet that is constrained not to move at the grid points. However, you're right that the linear solve will take a lot of time, especially since the matrices are not sparse. It's really up to you.




#5219073 Cubic interpolation on an irregular grid (triangle mesh)?

Posted by cadjunkie on 25 March 2015 - 09:06 AM

Yes you can, in a similar fashion to how you're doing it there. You'll need to have more data points along the edges and in the middle (10 points for a cubic interpolation), similar to a Bezier triangle:




It's going to be tricky since it'll be hard to guarantee a 10-point triangle. There might be some sort of adaptive method, but I'm not sure how to make it work with the mesh you've shown here.

#5217872 Explaining to an idiot (me) about mathematics concerning ocean waves

Posted by cadjunkie on 20 March 2015 - 09:37 AM

From reading Tessendorf's paper and all the other code examples out there (and doing an implementation in Matlab), I think I can try to offer some explanations.


It seems you're getting the correct \(H_0\) spectrum from the texture. My image (in grayscale) and images in other papers look similar to yours:




To get the height map, you need to use your texture values in the following formula and you can plug in the time parameter there as well:

\[H(k,t) = H_0(k,t)e^{i \omega t} + H_0^* e^{-i \omega t},\, \omega = g \left | k \right |\]

where \(H_0^*\) is the complex conjugate of \(H_0\). I think since your Y coordinate in the color vec3 seems to be the "imaginary" component of the complex number, then you can get the complex conjugate by just multiplying Y by -1. This \(H(k,t)\) is your height map. Here you can plug in your time \(t\) and get the animated displacements (at least I think it works; I didn't do the full animation). The decomposition into the spatial directions is given by:

\[D_x(k)=i\frac{k_x}{\left | k \right |}H(k),\,D_y(k)=i\frac{k_y}{\left | k \right |}H(k)\]


From there, you take the inverse FFT to get your spatial coordinates. In my Matlab implementation, those came out as complex numbers so I just took the norm and used that for the values, but I think you can just use the real parts as well. My 3 spectrum images look like this (dx, height, dy):




I haven't gotten around to using these to actually create waves. I've never done water simulations before, but the problem interested me from a theoretical standpoint so I figured I'd try it out. I'll post my Matlab code and if you use NumPy/SciPy the functions I used are very similar if you want to load it up and mess with it. Hope this helps!

close all;

M = 512;    % resolution of grid
N = 512;    % resolution of grid
L = 2000;   % height of waves used to calculate dx and dz

[m,n] = meshgrid(-M/2:1:M/2, -N/2:1:N/2);
kx = 2*pi/L * m;
ky = 2*pi/L * n;

W = [30,0];        % horizontal wind vector (m/s)
V = norm(W);        % m/s?
What = W / V;
g = 9.81;           % m/s^2
A = 1.0;            % amplitude of wave

% compute Phillips spectrum and freq amplitudes
P = zeros(M,N);
H0 = zeros(M,N);
for i=1:M
    for j=1:N
        k = [kx(i,j),ky(i,j)];
        nk = norm(k);
        if(nk == 0)
            mag = 0;
            khat = k / nk;
            mag = A / nk^4 * exp(-1/(nk * V^2/g)^2) * ((dot(khat,What)))^2;
        damping = exp(-nk^2*(L/1000)^2);
        val = mag * damping;
        P(i,j) = val; 
        H0(i,j) = 1/sqrt(2) * (randn + randn * 1i) * sqrt(val);


% compute heightmap
H = zeros(M,N);
Dx = zeros(M,N);
Dy = zeros(M,N);
t = 1;
for i=1:M
    for j=1:N
        k = [kx(i,j),ky(i,j)];
        nk = norm(k);
        w = sqrt(g*nk);
        val = H0(i,j)*exp(1i*w*t) + conj(H0(M-i+1, N-j+1))*exp(-1i*w*t);
        H(i,j) = val;
        if(nk == 0)
            Dx(i,j) = 0;
            Dy(i,j) = 0;
            Dx(i,j) = 1i*k(1)/nk * val;
            Dy(i,j) = 1i*k(2)/nk * val;

dx = abs(ifft2(Dx));
dy = abs(ifft2(Dy));
dz = abs(ifft2(H));


#5214574 How to do 3D feature detection?

Posted by cadjunkie on 04 March 2015 - 06:01 PM

I work directly with Parasolid and if it has feature recognition capabilities, I've never seen them. I know that other CAD packages such as Solidworks have tools to do that, but those tools may only exist within the CAD program and not within the geometry kernel itself.


This is actually a big area of research in the CAD world because there are standards such as IGES and STEP that just give B-spline geometry but no "feature" information. We call them "dumb bodies" because they don't tell you anything about how they were created. However, B-spline geometry can give you clues as to what certain shapes are. For example, a NURBS curve might be a circular edge based on how the control points relate to each other. But in this case, you have a curve you can test and not a string of vertices.


I'm sure that simple features such as holes, bosses, etc. can be extracted from the point data but it would be somewhat more involved. For example, the triangle vertices on a circular edge will have some roundoff error so you'll have to tolerance that, etc. I don't know offhand of any software tools that do this kind of thing that don't cost an arm and a leg.

#5213989 Delaunay Traingulation

Posted by cadjunkie on 02 March 2015 - 12:21 PM

I think you're better off using Fortune's algorithm. Not only do you get exactly what you're after, the time complexity is pretty good.

#5211497 How to create code to map function for Coefficient of Lift

Posted by cadjunkie on 18 February 2015 - 11:47 AM

My aircraft textbook states that the continuous function can be approximated by blending the functions between the coefficient of lift for a flat plate and a linear function that increases with angle of attack:


$$C_L = (1-\sigma(\alpha)) C_{L,linear} + \sigma(\alpha) C_{L,flat}$$




The flat plate function can be expressed as: $$C_{L,flat}=2\,\text{sign}(\alpha)\sin^2{\alpha}\cos{\alpha}$$


The linear function can be expressed as: $$C_{L,linear}=C_{L_0}+C_{L_{\alpha}}\alpha$$

$$C_{L_{\alpha}} = \frac{\pi AR}{1+\sqrt{1+(AR/2)^2}}$$



Here, b is the wingspan and S is the wing aspect ratio, and CL0 is a dimensionless quantity that depends on the aircraft but you can probably use 0 or something. The blending function you can use is:


$$\sigma(\alpha) = \frac{1+e^{-M(\alpha-\alpha_0)}+e^{M(\alpha+\alpha_0)}}{(1+e^{-M(\alpha-\alpha_0)})(1+e^{M(\alpha+\alpha_0)})}$$


Here, M is the transition rate and alpha_0 is the cutoff, so you can adjust the parameters to get the curve you want. When I was doing a small aircraft simulator for UAVs, I used these same equations so I know they work.