QuoteOriginally posted by: CuchulainnQuoteOriginally posted by: PolterQuoteOriginally posted by: CuchulainnQuoteFDM is an algorithm, there is no 'mathematical correctness' regarding element indices. Just convention and personal preference. Let the choice of prefered convention up to the user.Mathematicians start at index == 1, C at 0, Fortran at 1. (std::vector always at 0). FDM maths also starts at 1.That does not hold in general and it depends on convention, if anything.When writing an FDM solver for u(t,x) discretized as v(n,m), I was very explicitly treating t=0 as n=0 (the intuitive choice). The problem was interoping with software that couldn't do that (and the resulting mismatch between the notation and the given software -- not elegant). It's quite the same for a barrier at stock price x=0 (and thus m=0). If you have IC or BCs that include 0 in your domain, then it's the most natural choice.If anything, the situation gets different if you want domains with negative numbers (like solving a heat equation on x in [-1, +1]) or complex domains -- but then neither choice helps, IMHO. Then, one has to use an equivalent of std::begin and std::end for each axis -- and the internal representation/indexing of a container doesn't matter anyway.Can you give a simple 101 example of what you are describing, especially how begin() pans out in a fd scheme. It's not clear what the advantages would be. What if I want 'ghost points'??I have a Range class in both x and t which then get discretised. So, I distinguish between the PDE domain and the discrete FD domain. No problem. Suppose that the notation is as above -- u(t,x) is the solution function, v(n,m) is the discretization.Let's start with a situation where we fix the initial indexes at 0 (for you it'd be 1, doesn't matter here).Suppose the user has specified the domain of "u" by delimiters, for simplicity assume rectangular domain -- t_begin, t_end, x_begin, x_end.At some point, when filling the values given by IC and BCs we want to write at those subsets of the domain.Here, for simplicity assume we have a way to access the values function "u" takes there (IC and BCs are given) with an interface u(t_begin, x), u(t, x_begin), u(t, x_end) (an attempt to access an interior point here not given by IC or BCs is not applicable, might as well assume the behavior is undefined in such a case).At some point, our code will contain something like this:for all interior m: v(0,m) = u(t_begin, x); // ICfor all interior n: v(n,0) = u(t, x_begin); v(n,m_steps) = u(t, x_end); // BCsNow, t_begin in "u" might have been even negative, but it doesn't matter, since we're writing to container "v" anyway.Further, we can abstract from the initial indexes by introducing n_begin, n_end, m_begin, m_end (m_steps above now gives the # of steps from m_begin to m_end):for all interior m: v(n_begin,m) = u(t_begin, x); // ICfor all interior n: v(n,m_begin) = u(t, x_begin); v(n,m_end) = u(t, x_end); // BCsNow, one can probably notice that we have notions of begin/end for each axis -- both (t,x) and (n,m).As you have noticed, you introduce a notion of the range at this point.Now, as a concept, range supports a begin/end access interface to provide delimiting values (modulo dereferencing).So, if you have a discretized range for "m" ((here: first and only, but it doesn't matter) spatial axis of the discrete FD domain) being a discretization of the range for "x" ((here: first and only, but it doesn't matter) spatial axis of the PDE domain), you can use begin(m), end(x), etc... instead of the above.When I was writing my code, I couldn't do that in a standard way (no C++11 then, member functions won't work with everything), but now that we have a generic standard notion of those delimiters, perhaps we can support this (providing suitable overloads for non-C++11 compilers, of course) -- makes sense?
Last edited by Polter
on November 8th, 2011, 11:00 pm, edited 1 time in total.