Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I have the following finite difference operator: $$ Lu_{ijk}:= du_{ijk} +c(u_{i-1,j,k} + u_{i+1,j,k} + u_{i,j-1,k} + u_{i,j+1,k} + u_{i,j,k-1} + u_{i,j,k+1})\\ -u_{i-1,j+1,k}-u_{i-1,j-1,k} - u_{i-1,j,k-1} - u_{i-1,j,k+1} \\ -u_{i+1,j+1,k}-u_{i+1,j-1,k} - u_{i+1,j,k-1} - u_{i+1,j,k+1} \\ - u_{i,j+1,k-1} - u_{i,j-1,k-1} - u_{i, j+1,k+1} - u_{i,j-1,k+1}, $$ where $d$ and $c$ are known constants. This is way too complex of an expression to expect people to understand, so I want to represent it by a stencil:

5-point stencil

For example, I'd put a sphere of one color (representing $du_{ijk}$) at the origin, a sphere of another color at $(-1,0,0)$, (representing $cu_{i-1,j,k}$), and so on, and finally spheres of a final color at $(-1,1,0)$ (representing $-u_{i-1,j+1,k}$). My preliminary Mathematica commands are:

Graphics3D[{Specularity[White, 10], Red, Sphere[{0, 0, 0}],
            Blue, Sphere[{0, 0, 1}], Sphere[{1, 0, 0}], Sphere[{-1, 0, 0}],
            Sphere[{0, -1, 0}], Sphere[{0, 0, -1}],
            Green, Sphere[{1, 1, 0}], Sphere[{1, 0, 1}]}, Boxed -> False]

but these spheres are way too big and cover each other up. How can this be fixed?

share|improve this question
add comment

3 Answers 3

up vote 10 down vote accepted

Take it in steps:

  1. Extract the coefficients and locations into an appropriate data structure.

  2. Use that data structure to create the graphics.

By examining the FullForm of the original expression, we can cobble a rule to find the key data: the coefficients $c$, $d$, and $-1$ and the offsets to the indexes. First, the expression itself:

s = Subscript;
exp = d s[u, i, j, k] + 
  c (s[u, i - 1, j, k] + s[u, i + 1, j, k] + s[u, i, j - 1, k] + 
     s[u, i, j + 1, k] + s[u, i, j, k - 1] + s[u, i, j, k + 1]) - 
  s[u, i - 1, j + 1, k] - s[u, i - 1, j - 1, k] - s[u, i - 1, j, k - 1] - s[u, i - 1, j, k + 1] - 
  s[u, i + 1, j + 1, k] - s[u, i + 1, j - 1, k] - s[u, i + 1, j, k - 1] - s[u, i + 1, j, k + 1] - 
  s[u, i, j + 1, k - 1] - s[u, i, j - 1, k - 1] - 
  s[u, i, j + 1, k + 1] - s[u, i, j - 1, k + 1]

Expand pairs the coefficients with the subscripts and Cases extracts the essential information:

data = Cases[Expand@exp, Times[c_, s[u, i_, j_, k_]] :> {c, i, j, k}];
data // MatrixForm

$$\left( \begin{array}{cccc} -1 & -1+i & -1+j & k \\ -1 & -1+i & j & -1+k \\ c & -1+i & j & k \\ \cdots \\ -1 & 1+i & 1+j & k \end{array} \right)$$

We are really interested in the offsets to the central index $(i,j,k)$, so one more step to extract them (via Replace this time) will be helpful. After doing it, let's group the offsets by common coefficient using GatherBy:

spec = GatherBy[{First@#, 
         Replace[Rest@#,  {Plus[x_?NumberQ, i_] ->  x, x_Symbol -> 0}, 1]} & /@ data, First] 

To illustrate, here is what the first few elements of the first entry in spec look like:

$$\left( \begin{array}{cc} -1 & \{-1,-1,0\} \\ -1 & \{-1,0,-1\} \\ -1 & \{-1,0,1\} \\ \cdots \\ -1 & \{1,1,0\} \end{array} \right)$$

(You might be happier just entering the data in this format, or something close to it, at the outset: it's easier than entering all those subscripts.)

Choose some colors:

colors = Array[Hue[# / Length@spec, .8, .8] &, Length@spec];

The rest is easy. Let's make sure to include some visual cues such as thin lines connecting the base point to its neighbors, for otherwise this will look only like a random jumble of balls.

Graphics3D[ { 
  Table[{Specularity[White, 10], 
    GrayLevel[0.7], Tube[{{0, 0, 0}, Last@#}, 0.025] & /@ spec[[i]],
    colors[[i]], Sphere[Last@#, .2] & /@ spec[[i]]}, {i, 1, Length@spec}]}, 
 Boxed -> False, Axes -> True, AxesLabel -> {"i", "j", "k"} ]

(I leave the creation of a color key as an exercise :-).)

Figure

share|improve this answer
    
If I can nitpick... I prefer horizontal and vertical connections so that it is clearer which point is located where w.r.t. the center. –  rubenvb Jul 10 at 8:30
    
@rubenvb That's a great idea. How would you implement it? In this diagram, renderings of such connections would all run into one another. I can think of some possibilities, such as introducing a light 3D grid, or a set of nested transparent cubes, or other such visual references, but perhaps you have a specific form of visualization in mind? –  whuber Jul 10 at 13:10
    
no idea really, I am terrible at drawing with Mathematica. I don't see how they would run into each other though. Something like this comes to mind though. –  rubenvb Jul 10 at 13:56
    
@rubenvb Thank you: I see what you mean. There are arbitrary elements to that diagram--it could be drawn in many equivalent ways. However, something like it could be implemented in Mathematica, provided we settle on rules to determine which sets of line segments to draw. –  whuber Jul 10 at 14:06
add comment

Although you haven't exactly asked this, you might like to generate your graphic automatically by applying pattern matching on your difference operator. The basic idea is as below:

Clear[i, j, k];

op = Plus @@ 
  MapThread[Subscript[u, i - #1, j - #2, k - #3] &, 
   RotateRight[{0, -1, 1, 0, 0, 0, 0}, #] & /@ {0, 2, 4}]

Giving $op = u_{i-1,j,k}+u_{i,j-1,k}+u_{i,j,k-1}+u_{i,j,k}+u_{i,j,k+1}+u_{i,j+1,k}+u_{i+1,j,k}$

Graphics3D[
 Level[op /. 
     Subscript[_, i_, j_, k_] :> Sphere[{i, j, k}, .1] /. {i -> 0, 
     j -> 0, k -> 0}, {-3}] /. {s : Sphere[{0, 0, 0}, _] :> 
    Sequence[Red, s], s_Sphere :> Sequence[Blue, s]}, Axes -> True]

3D Stencil?

Disclaimer: I'm a fairly new and rather unsophisticated Mathematica user, so my code is probably crap and likely to break if you breathe on it too hard, but anyway, there you have it.

share|improve this answer
add comment

Not sure if this is what you're after. Anyway:

r = .1; 
Graphics3D[{Specularity[White, 10], 
            Red, Sphere[{0, 0, 0}, r], 
            Blue, Sphere[{0, 0, 1}, r], Sphere[{1, 0, 0}, r], Sphere[{-1, 0, 0}, r],  
                                        Sphere[{0, -1, 0}, r], Sphere[{0, 0, -1}, r],
            Green, Sphere[{1, 1, 0}, r], Sphere[{1, 0, 1}, r]}, Boxed -> False, 
            Axes -> True]

enter image description here

share|improve this answer
    
Oh, that was really easy! Thanks! –  Nick Thompson May 1 '13 at 23:15
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.