Hot answers tagged

7

Like WolfgangBangerth, I strongly recommend that you reconsider your motivation (or your supervisor's motivation) for this goal. First, look at Center for Exascale Simulation of Combustion in Turbulence (ExaCT) web site. They have some proxy applications that may be of use to you in testing your mechanisms, and publications by authors such as John Bell, Joe ...


7

You have a size mismatch issue: A is a count x 2*count matrix, and you are trying to solve Ax=B with B a 2*count x 1 vector. Moreover, if you compile without -DNDEBUG, you should get a nice assertion telling you this is wrong. To resize a matrix or vector: A.resize(count, 2*count); B.resize(count);


7

The formula still applies for vectors - I find it helpful to use index notation. If $(\nabla u) = u_{,i}$, where $i$ implicitly refers to the $i$th entry of the gradient (i.e. $\frac{\partial u}{\partial x_i}$), and we assume Einstein notation (i.e. repeated indices sum over that index), then we can rewrite $\nabla \cdot \nabla {\bf u} = (u_j)_{,i,i} = \...


6

A cross product will tell you both the magnitude and the sign of your angular velocity. In general, angular velocity is defined by a vector as $$\vec{\omega} = \dfrac{\vec{r}\times\vec{v}}{|\vec{r}|^2}$$ where $\vec{\omega}$ is the angular velocity, $\vec{r}$ is the vector from the center of rotation to the point under consideration, and $\vec{v}$ is the ...


4

This is a O(m) problem where m is the size of the big array. So this should be fast: int n = v[v.size()-1]; int k=0; VectorXi lengths[n+1], starts[n+1]; for(int i=0; i<v.size(); ++i) { int i0 = i; while(v[i]==k && i<n) ++i; lengths[k] = i-i0; starts[k] = i; ++k; }


4

What you suggest is a multi-year research program. Don't even think about doing this yourself, you will not within several years achieve anything that comes close to the current research. Rather, search the web for codes that are available and adjust them to your needs. A good first place would be the Department of Energy's Computational Combustion research ...


4

Certainly. There are a few things you have to define for your type that are listed on this page in the documentation: https://eigen.tuxfamily.org/dox/TopicCustomizing_CustomScalar.html It basically boils down to defining arithmetic operators appropriately for your type, plus specializing a traits template NumTraits that describes your type. The link above ...


3

You can calculate the stream function yourself, and from it you can draw contours or streamlines of constant $\psi$. Let us assume you are on a two-dimensional incompressible flow $\mathbf{u}=(u,v,0)^T$, then you can find an exact differential $d\psi$ that satisfies the mass conservation equation: $ \frac{\partial u}{\partial x } + \frac{\partial v}{\...


3

My interpretation of this sentence is the same as yours, and I would have implemented the algorithm the same way.


3

The short answer is "It depends". Certainly, you can try to find the Helmholtz decomposition on your sampled data, and find your irrotational and solenoidal components. However, there are certain requirements on your original vector field you started with. In general, it is that the vector field you are trying to decompose has to be sufficiently smooth and ...


2

For plotting, it is easier in my opinion to not use meshgrid if you want to scale the arrows. You have a vector field $(E_X, E_Z)$ and you can simply normalize it like in the code below: clear; clc; p = 1; e = 8.85*10^(-12); x =linspace(-5 , 5, 20); z = linspace(-5 , 5 ,20); for i = 1:length(x) for k = 1:length(z) R=sqrt(x(i)^2 + z(k)^2) ; ...


2

you just need to add output to your min() function, e.g.: [value, ind] = min(v) ind is then the index values of your min locations (the i,j that you are looking for) and value has the corresponding min values


2

You're correct, if the orientational vector is unitary. Otherwise you must calculate the unitary vector in the direction of $\vec{O}_1$ and then perform the projection (just divide by the norm of the orientational vector). The perpendicular component is just the rejection vector of $\Delta \vec{R}$, given by the definition: $\Delta\vec{R}_{\perp} = \Delta\...


2

I will try to answer the part of the question regarding the accuracy of the calculation, which certainly affects all the other things. Integrals of the type ($T_j$ denoting the $j$th triangle in the discretization and $\vec{x}$ some arbitrary observation point) $$ \int\limits_{T_j}\frac{1}{||\vec{x}-\vec{y}||}dy,\quad \vec{y}\in T_j $$ have a singularity ...


2

You can think of each vector as a point in your linear space. As such, we can use a simple quadtree/octree-like algorithm to map your points into boxes, with "nearby" vectors assigned to the same or an adjacent box. With $n$ total vectors the vector-to-box map costs $\mathcal{O} (n\log n)$, and once this is done you can choose $m$ boxes and select the vector ...


2

Easy: Decompose the polygons into (unoriented) line segments each of which is sorted by vertex index: $$ A = \{[1,2], [2,3], [3,4], [4,5], [1,5]\}, \\ B = \{ [1,6], [6,7],[7,8],[3,8], [2,3], [1,2] \}. $$ Then you want to want to consider the union of all of these edges, but removing the duplicated ones. So you need to form $$ (A \cup B) \setminus (A \...


2

What you're describing is also a critical step in the k-nearest-neighbours method. So no need to reinvent the wheel, we can just look how other people have sped up that algorithm. I don't know about any dictionary like structure that returns this directly, but you could use a k-d tree. If properly implemented, you can get the k closest vectors pretty quickly....


2

Of course you can! The starting vector is completely arbitrary. You'll get breakdown at step 1 if you do though. Be sure to check what it entails. This will typically happen by a "happy accident" rather than by your explicit choice: if you already know an eigenvector then you don't need to run Lanczos to compute it (or you need to run it with a ...


1

For anyone interested in this problem, I found the following solution: #include "itensor/all.h" using namespace itensor; int main() { int N = 100; // // Initialize the site degrees of freedom. // auto sites = SpinHalf(N,{"ConserveQNs=",false}); //make a chain of N spin 1/2's //Transverse field Real h = 4.0; // //...


1

Ok I will answer my own question: The problem is in the line solve(L==a,A,bc) which needs to be replaced by solve(a==L,A,bc). The two versions seem to non interchangeable. Doing so will result in the magnetic field B=project(curl(A),V) which is shown in the picture below:


1

Somehow I have to use t-sne, but I really don't know how. Since you have a PhD my answer will be brief, it's quite a lengthy subject. The aim of dimensionality reduction is to preserve as much of the significant structure of the high-dimensional data as possible in the low-dimensional map. For high-dimensional data that lies on or near a low-dimensional, ...


1

If you know how the x-axis relates from one dataset to the next, the best solution I can think of would be to use the command 'interp1' in matlab. I've compiled a small script with output to show you some example usage. Essentially what you want to do is remap your data from each dataset so that they share some common axis, which you can do with 'interp1' in ...


1

Initialize the cluster centers as your subset $V\in X$, where $V=\{x_i\}$. Then run a couple of K-mediods iterations. After that you will see that the certain vectors will come closer, essentially trying to represent similar peaks. It is then possible to merge them. Another way to do this is mean-shift algorithm, where the modes are being updated (I think ...


1

Sounds like you want to thin your data where it is dense, and learn the support of your data summarized by data points. If you don't have too many points, you can generate a distance matrix, and prune the points with the closest neighbors. (I don't think this method has a name.) Otherwise, if you are trying to find the extremal points, then archetypal ...


1

After some searching I found a C/C++ function called project() in the article 'Real-Time Fluid Dynamics for Games' by Jos Stam http://www.dgp.toronto.edu/people/stam/reality/Research/pdf/GDC03.pdf I had to do a little bit of editing, such as increasing iteration variable k (dependent on how divergent your input field is) and employing periodic boundaries. ...


1

While I don't state this to be an optimal implementation or one taking into account the all edge cases, below is a sample implementation I put together for your task that can give you an idea of how one might approach your problem. The Header File #ifndef _container_hpp_ #define _container_hpp_ #include <vector> class container { public: class ...


1

From your comment "The plot emphasises expectedly the connections of node 1", I guess that maybe the idea is showing that node 1 is connected not only to a set of nodes that are "one close to the other", but its interconnections are spread somewhat evenly across the graph. Reverse Cuthill-McKee in theory reorders the nodes so that clusters are mapped into ...


1

If you just want an approximation to the sum of their norms, you can choose a random sample of $m$ of the indices (i.e., a random set $S \subseteq \{1,2,\dots,N\}$ of size $m$) and compute $$f(\textbf{x}) = \sum_{n \in S} \| x_n \|_2.$$ By a standard central limit theorem argument, $f(\textbf{x})$ will be a reasonable approximation to the sum of the norms; ...


1

Is it possible to tell in index notation whether a vector is a row or column vector, or is that supposed to be clear based on context? It seems the answer is actually lurking in your question itself here: if $u \cdot v \equiv u^T v$ then that doesn't leave much wiggle room. $u^T$ will have to be a row vector in order for the resulting product to be a ...


Only top voted, non community-wiki answers of a minimum length are eligible