Register a SA Forums Account here!
JOINING THE SA FORUMS WILL REMOVE THIS BIG AD, THE ANNOYING UNDERLINED ADS, AND STUPID INTERSTITIAL ADS!!!

You can: log in, read the tech support FAQ, or request your lost password. This dumb message (and those ads) will appear on every screen until you register! Get rid of this crap by registering your own SA Forums Account and joining roughly 150,000 Goons, for the one-time price of $9.95! We charge money because it costs us money per month for bills, and since we don't believe in showing ads to our users, we try to make the money back through forum registrations.
 
  • Post
  • Reply
nielsm
Jun 1, 2009



How to read your error message:

quote:

code:
Fit.c: In function ‘main’:
Fit.c:77:1: error: incompatible type for argument 1 of ‘Lorentzian’
Fit.c:9:6: note: expected ‘float’ but argument is of type ‘float *’
The lines all start with "Fit.c" so that's the file the error is occurring in. The first line then tells you that these things are happening while compiling the "main" function.
So it could help if you posted your entire main() function exactly as-is.

The following lines begin with "Fit.c:XX:YY:", that should be read as "line XX, position YY in Fit.c". So it's telling you that lines 77 and 9 are related to the error.
Taking the function prototype...
C++ code:
void Lorentzian(float phi, float *alpha, float *nu, float *nunought, int index)
The first parameter is declared to be of type "float", but according to line 3 in the error message, you are passing a float* (pointer to float). You need to not pass a pointer, but I can't tell you what you should be doing instead since you haven't posted your call site.

The three next errors (Note that each error takes up two lines in the message: one telling you about the type of error, and one informative telling you the formal and actual types involved.) all mention expecting a float* but getting a float*[size], presumably here you are passing an array directly. It can work, but you are better off being explicit about it. Maybe &nu[0] instead of just nu. That gives you a pointer to the first element of the array.

Take a shot yourself at the final two errors for that call.

Adbot
ADBOT LOVES YOU

rjmccall
Sep 7, 2007

no worries friend
Fun Shoe
E;FB

Congratulations on quoting a ton of code without managing to show us the lines with all the errors on them, i.e. the places where you're calling Lorentzian and fit.

code:
Fit.c: In function ‘main’:
Fit.c:77:1: error: incompatible type for argument 1 of ‘Lorentzian’
Fit.c:9:6: note: expected ‘float’ but argument is of type ‘float *’
You're probably passing phi as an argument instead of some element of it.

code:
Fit.c:77:1: warning: passing argument 2 of ‘Lorentzian’ from incompatible pointer type [enabled by default]
Fit.c:9:6: note: expected ‘float *’ but argument is of type ‘float (*)[(sizetype)(size)]’
Fit.c:77:1: warning: passing argument 3 of ‘Lorentzian’ from incompatible pointer type [enabled by default]
Fit.c:9:6: note: expected ‘float *’ but argument is of type ‘float (*)[(sizetype)(size)]’
Fit.c:77:1: warning: passing argument 4 of ‘Lorentzian’ from incompatible pointer type [enabled by default]
Fit.c:9:6: note: expected ‘float *’ but argument is of type ‘float (*)[(sizetype)(size)]’
For each of these, you're probably passing &alpha, where that's some variable-length local array. You should just pass alpha instead.

code:
Fit.c:77:1: warning: passing argument 6 of ‘fit’ from incompatible pointer type [enabled by default]
Fit.c:21:6: note: expected ‘int *’ but argument is of type ‘float *’
This is probably just you passing the wrong thing.

code:
Fit.c:77:1: warning: passing argument 9 of ‘fit’ from incompatible pointer type [enabled by default]
Fit.c:21:6: note: expected ‘char *’ but argument is of type ‘char **’
This is probably you passing the address of something instead of that thing.

Shugojin
Sep 6, 2007

THE TAIL THAT BURNS TWICE AS BRIGHT...


Oh god drat it I did forget the call :downs:

code:
fit(nu, phi, sigma, n, a, ia, ma, *Lorentzian(phi, &nu, &a, &nunought, n),&dunnowhatthisdoes);
where &dunnowhatthisdoes is defined as

code:
char* dunnowhatthisdoes;

nielsm
Jun 1, 2009



Right. Looking at the declaration of fit() again:
C++ code:
void fit(float nu[],float phi[],float sig[],int n,float a[],int ia[],int ma,
          void (*func)(float,float *,float *,float *,int),char *label)
It takes 3 float arrays, an int n, two more arrays (a float one and an int one), yet another int, then a function pointer, and finally a C string.

Your primary mistake is that you are calling Lorentzian() instead of passing it as a pointer. The fit() function itself will be the one calling it (repeatedly, presumably) to do its work, so you shouldn't call it yourself, particularly not at the call site for fit().

Try this:
C++ code:
float a[other_size];
int ia[yet_another_size];
int ma = 42;
char *label = "whatever";
fit(nu, phi, sigma, size, a, ia, ma, &Lorentzian, label);
Go read the documentation for that function so you actually know what it wants for the a and ia arrays, their purpose and what sizes they should have, as well as what the purpose of the label thing is.

tractor fanatic
Sep 9, 2005

Pillbug

GrumpyDoctor posted:

Assuming that this is workable, what's the most elegant way to account for the fact that these will actually be consumed by C++? (I'm not exposing C++ so that clients can use platforms other than the one I happened to compile with.) The way to define structures is different, so as I see it I can either have two parallel header files (ugh) or some kind of macro gadgetry controlled by #ifdef __cplusplus that will expand into either the C definition or the C++ one (also ugh). Is there some simpler way I'm missing?

#ifdef __cplusplus
extern "C"{

? I don't think C++ will complain about the typedef struct syntax.

ExcessBLarg!
Sep 1, 2001

GrumpyDoctor posted:

How workable is it to expose a C API that includes structures?
It's fine for the API, if your intention is for all the relevant components to be compiled from source using the same compiler (i.e., no binary plugins). Otherwise you're going to run into issues with the ABI. Here's some tips:
  • Define the struct in C, use 'extern "C"'.
  • Use stdint.h width-specified integer types (e.g., int32_t) instead of "int" or "long".
  • Use GCC's "aligned" or "packed" type attributes to enforce a known structure packing. If you want to retain compatibility with MSVC, use you can use structure-packing pragmas instead. Unfortunately there's no standards compliant way to do this.
  • Note that certain platforms (e.g., SPARC) can't access unaligned data elements and OS traps to emulate them are slow. You may want to use your ABI struct purely for marshaling.
You could also consider using an existing data serialization format, like XDR, if you want to ensure something that's portable across platforms and compiler versions without worrying over the gory details. It's also useful if the API representation of your data is ever to be stored in a file and possibly used cross-platform, and/or it may grow into a network/RPC-based API.

rjmccall
Sep 7, 2007

no worries friend
Fun Shoe
You really shouldn't need to screw around with portable-size typedefs and structure packing just to get interoperation between compilers on the C level; that's all supposed to be dictated by the platform ABI, and if compilers don't agree on that, they're not really targeting the same platform, and you can't really assume *any* level of interoperation. That is, if two compilers don't agree on how big an int is, there's not really any reason to assume that they'll pass one to a function the same way.

It gets weird on Windows, of course, because everything does, but elsewhere, there's really no reason to anticipate problems passing structs between, say, suncc and GCC or GCC and clang.

Now, if you're talking about using a struct as a wire format, e.g. for a binary file or a network protocol, then yes, you need to pull out all the portable layout stops.

raminasi
Jan 25, 2005

a last drink with no ice

tractor fanatic posted:

#ifdef __cplusplus
extern "C"{

? I don't think C++ will complain about the typedef struct syntax.

You know, I just kind of assumed that it wouldn't work but you're right, C++ doesn't give a poo poo about the "C" syntax :cripes:


rjmccall posted:

You really shouldn't need to screw around with portable-size typedefs and structure packing just to get interoperation between compilers on the C level; that's all supposed to be dictated by the platform ABI, and if compilers don't agree on that, they're not really targeting the same platform, and you can't really assume *any* level of interoperation. That is, if two compilers don't agree on how big an int is, there's not really any reason to assume that they'll pass one to a function the same way.

It gets weird on Windows, of course, because everything does, but elsewhere, there's really no reason to anticipate problems passing structs between, say, suncc and GCC or GCC and clang.

Now, if you're talking about using a struct as a wire format, e.g. for a binary file or a network protocol, then yes, you need to pull out all the portable layout stops.

This kind of thing is exactly what I was looking for. None of this will be used on the wire or serialized in any way, so I guess that cuts down on some of my problems. As it turns out, the primary sort of interop I'm expecting is actually MSVC <-> Other Versions of MSVC, which sounds like I should probably be okay. I vaguely recall that Microsoft changed stuff for 2010 but I'm trying to target a level that's lower than that. (I mean, "C compiled with 2010 can't talk to C compiled with 2008" sounds colossally stupid, even for Microsoft.)

Honestly, if my strategy has to be "assume it works in some cases and assume it doesn't work in others," that's fine, because I don't really have any other options. I just want to do what I can to minimize problems.

Plorkyeran
Mar 22, 2007

To Escape The Shackles Of The Old Forums, We Must Reject The Tribal Negativity He Endorsed
The only issue with mixing code compiled with multiple versions of visual studio is that it can lead to multiple C runtimes in your application, which can be exciting.

raminasi
Jan 25, 2005

a last drink with no ice

Plorkyeran posted:

The only issue with mixing code compiled with multiple versions of visual studio is that it can lead to multiple C runtimes in your application, which can be exciting.

I can't think of a way around that, so not my problem :pseudo: But seriously, I'm providing my own allocating and freeing functions, and these structures are all essentially big bundles of chars and doubles so I don't anticipate problems on that front.

MrMoo
Sep 14, 2000

Plorkyeran posted:

The only issue with mixing code compiled with multiple versions of visual studio is that it can lead to multiple C runtimes in your application, which can be exciting.

It would appear to be popular with many games, Steam installs love to pack various CRTs from 2005 upwards. Example being Borderlands which uses an old audio library.

FamDav
Mar 29, 2008

GrumpyDoctor posted:

You know, I just kind of assumed that it wouldn't work but you're right, C++ doesn't give a poo poo about the "C" syntax :cripes:

That's because the C syntax is perfectly valid.

ex.

C++ code:
typedef struct S_t {} S;

int main()
{
  S x;
  S_t y;
  struct S_t z;
}
is perfectly valid. Just realize that C++ implicitly aliases struct S_t by S_t.

What's interesting is that

code:
typedef struct S {};

typedef struct T {} T;
typedef T T;

int main()
{
  S x;
  T y;
}
Is also syntactically correct, which comes directly from the standard saying that the alias is optional in typedefs and that redefinition the same alias is allowed. You'll probably get a warning about having a typedef that declares no alias, though.

EDIT:

Oh, but the redefinition stuff does not work within classes, so you can't write something crazy like:

code:
template <typename T>
class assert_alias
{
  typedef T type;
  typedef type type;
};

FamDav fucked around with this message at 04:40 on Mar 7, 2013

Bruegels Fuckbooks
Sep 14, 2004

Now, listen - I know the two of you are very different from each other in a lot of ways, but you have to understand that as far as Grandpa's concerned, you're both pieces of shit! Yeah. I can prove it mathematically.

MrMoo posted:

It would appear to be popular with many games, Steam installs love to pack various CRTs from 2005 upwards. Example being Borderlands which uses an old audio library.

There was a time when I had to use three different versions of visual studio to get different dlls in a project I was working on to compile.

Oh yeah, it also had leaks like this one:
http://support.microsoft.com/kb/555563

Have fun.

raminasi
Jan 25, 2005

a last drink with no ice

FamDav posted:

That's because the C syntax is perfectly valid.

That's what I meant with the scare quotes.

Hammerite
Mar 9, 2007

And you don't remember what I said here, either, but it was pompous and stupid.
Jade Ear Joe
This is a header file from something I'm working on.

code:
#ifndef COORDINATE_RANGE_H
    #define COORDINATE_RANGE_H

    #include <utility>

    class coordinateRange {
        protected:
            int minimum;
            int maximum;
        public:
            coordinateRange (int const rMin, int const rMax);
            static coordinateRange checked (int limit1, int limit2);
            bool inRange (int const value) const;
            int length () const;
            bool isEmpty () const;
            int getMin () const;
            int getMax () const;
            void setMin (int const newMin);
            void setMax (int const newMax);
            void setMinChecked (int const newMin);
            void setMaxChecked (int const newMax);
            void makeEmpty ();
            coordinateRange intersection (coordinateRange const that) const;
            coordinateRange unionHull (coordinateRange const that) const;
            void restrictXforYbyZ (int const y, coordinateRange const zRange);
            void restrictXforYbyA (int const y, coordinateRange const aRange);
            void restrictXforYbyB (int const y, coordinateRange const bRange);
            void restrictXforYbyC (int const y, coordinateRange const cRange);
            void restrictXforZbyY (int const z, coordinateRange const yRange);
            void restrictXforZbyA (int const z, coordinateRange const aRange);
            void restrictXforZbyB (int const z, coordinateRange const bRange);
            void restrictXforZbyC (int const z, coordinateRange const cRange);
            void restrictYforXbyZ (int const x, coordinateRange const zRange);
            void restrictYforXbyA (int const x, coordinateRange const aRange);
            void restrictYforXbyB (int const x, coordinateRange const bRange);
            void restrictYforXbyC (int const x, coordinateRange const cRange);
            void restrictYforZbyX (int const z, coordinateRange const xRange);
            void restrictYforZbyA (int const z, coordinateRange const aRange);
            void restrictYforZbyB (int const z, coordinateRange const bRange);
            void restrictYforZbyC (int const z, coordinateRange const cRange);
            void restrictZforXbyY (int const x, coordinateRange const yRange);
            void restrictZforXbyA (int const x, coordinateRange const aRange);
            void restrictZforXbyB (int const x, coordinateRange const bRange);
            void restrictZforXbyC (int const x, coordinateRange const cRange);
            void restrictZforYbyX (int const y, coordinateRange const xRange);
            void restrictZforYbyB (int const y, coordinateRange const aRange);
            void restrictZforYbyC (int const y, coordinateRange const bRange);
            void restrictZforYbyA (int const y, coordinateRange const cRange);
    };
#endif
The last 24 methods are to restrict the range to its intersection with a supplied range in another coordinate (for a fixed value of a third coordinate). I feel like I'm missing some conceptual tool that would allow me to write fewer methods, however I don't know enough to know what it is. I worry I'm probably writing code that will be ill-performing as well as bloated with methods.

Perhaps I should have a coordinateRange carry around information about which coordinate it is "in". They could be subclasses, and I could implement virtual methods for the restriction operation. But I don't know if that's making things too complicated.

rjmccall
Sep 7, 2007

no worries friend
Fun Shoe
Well, one possibility would be to make distinct X, Y, and Z coordinate types, make coordinateRange be generic over them, and then just use overload resolution.

Spiritus Nox
Sep 2, 2011

Hey there. I'm working on an algorithm that reads in a pair of integers x and y (assuming x and y are relatively prime and part of the Stern-Brocot Tree - x and y are also handled as doubles in order to more accurately determine x/y), and reports the 'turns' one would have to take along the tree to find x/y - for example, passing in 5 and 7 returns "LRRL" - a left turn, two right turns, and another left turn on the Stern-Brocot Tree. The upload site that grades my program reports an output error (IE my program is producing a wrong answer), but I've yet to find any test data that breaks my algorithm. If anyone has any insight, I'd be very grateful. Here's the algorithm:
code:

string findSternBrocot(double x, double y)
{
    string s = "";
    double num = 1, denom = 1;
    double current = num/denom;
    double target = x/y;
    double leftNum = 0, leftDenom = 1;
    double rightNum = 1, rightDenom = 0;
    bool done = false;

    while(!done)
    {
        if(current < target)
        {
            leftNum = num;
            leftDenom = denom;
            num += rightNum;
            denom += rightDenom;

            current = num/denom;
            s += "R";
        }
        else if(current > target)
        {
            rightNum = num;
            rightDenom = denom;

            num += leftNum;
            denom += leftDenom;

            current = num/denom;
            s += "L";
        }
        else
        {
            if(num == denom)
            {
                s += "I";
            }
            done = true;
        }
    }

    return s;
}

Spiritus Nox fucked around with this message at 22:27 on Mar 7, 2013

tractor fanatic
Sep 9, 2005

Pillbug
I don't see anything wrong with your algorithm, although I'd use pure integer math for it. What about the rest of your program?

Jabor
Jul 16, 2010

#1 Loser at SpaceChem
E: n/m, misread the code.

Have you tried your solution on very large pairs of numbers?

Jabor fucked around with this message at 23:04 on Mar 7, 2013

Spiritus Nox
Sep 2, 2011

tractor fanatic posted:

I don't see anything wrong with your algorithm, although I'd use pure integer math for it. What about the rest of your program?

Not much to tell, but here it is:

code:
int main()
{
    //Data Abstraction
    double l, r;
    bool done = false;
    string s;

    while(!done)
    {
        //Input
        cin >> l >> r;

        if(l == 1 && r == 1)
        {
            done = true;
        }
        else
        {
            //Process: The 'turns' required to traverse the tree and find the
            //         supplied stern brocot number are calculated and returned
            //         by the findSternBrocot function, discussed in greater detail
            //         below.
            s = findSternBrocot(l, r);

            //Output
            cout << s << endl;
        }


    }
Just reads in two numbers, then outputs the sequence of turns - again, the assignment assumes that l/r is contained within the Stern Brocot Tree. 1/1 is our terminator case, so it's not processed.

Also, I'm using doubles because, as you can see, I'm determining my 'target' point in the tree by finding the decimal value of l/r, so I figured it would be better not to lose precision like I would with integers.


Jabor posted:

Have you tried your solution on very large pairs of numbers?

Largest input I've tested is 200000/200003. I've tried larger, adding too many more digits starts producing such massive walls of Ls and Rs that it becomes really hard to tell if my output differs from that of the upload site.

Nippashish
Nov 2, 2005

Let me see you dance!

Spiritus Nox posted:

Also, I'm using doubles because, as you can see, I'm determining my 'target' point in the tree by finding the decimal value of l/r, so I figured it would be better not to lose precision like I would with integers.

Checking a/b < c/d is equivalent to checking a*d < c*b and the second can be done with only integers.

Spiritus Nox
Sep 2, 2011

All right, I've tried implementing it that way and I'm currently getting an infinite loop - I'm probably missing something really basic, but here's my code thus far.

code:
string findSternBrocot(int l, int r)
{
    string s = "";
    int targetNum = l, targetDenom = r;
    int num = 1, denom = 1;
    int current = num*targetDenom;
    int target = targetNum*denom;
    int leftNum = 0, leftDenom = 1;
    int rightNum = 1, rightDenom = 0;
    bool done = false;

    while(!done)
    {
        if(current < target)
        {
            leftNum = num;
            leftDenom = denom;
            num += rightNum;
            denom += rightDenom;

            current = num*targetDenom;
            s += "R";
        }
        else if(current > target)
        {
            rightNum = num;
            rightDenom = denom;

            num += leftNum;
            denom += leftDenom;

            current = num*targetDenom;
            s += "L";
        }
        else
        {
            if(num == denom)
            {
                s += "I";
            }
            done = true;
        }

        cout << num << "/" << denom << " " << current << " " << target << " " << s << endl; //Debug purposes
    }

    return s;
}

OneEightHundred
Feb 28, 2008

Soon, we will be unstoppable!
I'm not sure what you're passing in there, but make sure that you're not overflowing beyond the integer resolution. You mentioned 200000/200003, multiplying those gives you 40000600000, but anything over 134217727 is an overflow and the higher bits will be truncated.

Spiritus Nox
Sep 2, 2011

OneEightHundred posted:

I'm not sure what you're passing in there, but make sure that you're not overflowing beyond the integer resolution. You mentioned 200000/200003, multiplying those gives you 40000600000, but anything over 134217727 is an overflow and the higher bits will be truncated.

I'm getting an infinite loop on 5/7, so overflow shouldn't be the problem there.

Hammerite
Mar 9, 2007

And you don't remember what I said here, either, but it was pompous and stupid.
Jade Ear Joe

rjmccall posted:

Well, one possibility would be to make distinct X, Y, and Z coordinate types, make coordinateRange be generic over them, and then just use overload resolution.

So... 6 classes xCoord, yCoord, zCoord, aCoord, bCoord, cCoord

A template coordinateRange, used like coordinateRange<yCoord>

Perhaps a template like lineSegment (a coordinateRange through a point)

Isn't there an alternative, to make an abstract "coordinate" class and have the 6 other classes implementations of it? Then things like coordinateRange would extend the other classes I guess.

Spiritus Nox
Sep 2, 2011

I got my program working - I just needed to update Target as well as Current in my modified algorithm. It passes now. Thanks for the tips!

Shugojin
Sep 6, 2007

THE TAIL THAT BURNS TWICE AS BRIGHT...


HILARIOUSLY dumb problem, I hope!

I'm trying to determine what the poo poo I want to do with my Lorentzian (and/or Gaussian) functions and I've encountered something annoying.

Here's a sample Lorentzian:
code:
void Lorentzian(float nu, float *a, float *phi, float *dyda, int n)
{
*phi = (1/PI)*(a[0]/(pow((nu - a[1]),2) + (pow(a[0], 2))));
dyda[0] = ((pow((a[1]-nu),2) - pow(a[0],2))/(PI*pow(pow(a[0],2)*(pow((a[1] - nu),2)),2))); /* derivative wrt alpha, a[0] = alpha, a[1] = nu 0*/
dyda[1] = (2*a[0]*(nu - a[1]))/(pow((pow(a[0],2) + pow((a[1] - nu),2)),2)); /*derivative vrt nu 0 a[0] = alpha, a[1] = nu 0*/ 

}
Heeeere's a little for loop!

code:
for (i=0;i<=size;i++)
	{
	Lorentzian(nu[i],a,phi,dyda,n);
	printf("%f\t %f\t %f\n",phi[i], dyda[0],dyda[1]);
	}
Now that works and prints poo poo to a screen but I really want it in a little debug file because I find that more convenient (because I can save older debug files for reference much more easily). That isn't working. Here's the relevant excerpts:

code:
FILE *debug;
debug = fopen("debug.txt", "w+");

for (i=0;i<=size;i++)
	{
	Lorentzian(nu[i],a,phi,dyda,n);
	fprintf(debug,"%f\t %f\t %f\n",phi[i], dyda[0],dyda[1]);
	}
WHY DOESN'T THIS PRINT TO A FILE AT ALL AUGH I CANNOT SEE IT

tractor fanatic
Sep 9, 2005

Pillbug
Close the file.

Shugojin
Sep 6, 2007

THE TAIL THAT BURNS TWICE AS BRIGHT...


:doh: thank you.

rjmccall
Sep 7, 2007

no worries friend
Fun Shoe

Hammerite posted:

So... 6 classes xCoord, yCoord, zCoord, aCoord, bCoord, cCoord

Side note: I have no idea what you're trying to model here, and the signatures of these restrict methods make no sense to me. Is this a six-dimensional space? How do you restrict a coordinate range by a range in a different dimension, and what does the fixed point in the third dimension have to do with anything? It's hard, and probably fundamentally misguided of me, to suggest abstractions when I don't understand the problem space, but I'll keep trying anyway.

Anyway, you would obviously want these classes to share an implementation, instead of just copying them N times; among other things, that makes it a lot easier to abstract over them. I'd suggest:

code:
enum Axis {X,Y,Z,A,B,C};
template <Axis> class Coord {
  ...
};
typedef XCoord Coord<X>;
...

Hammerite posted:

A template coordinateRange, used like coordinateRange<yCoord>

Either that or coordinateRange<Y>.

Hammerite posted:

Perhaps a template like lineSegment (a coordinateRange through a point)

Sure.

Hammerite posted:

Isn't there an alternative, to make an abstract "coordinate" class and have the 6 other classes implementations of it? Then things like coordinateRange would extend the other classes I guess.

This is certainly an alternative. I don't know if it makes solving your problems any easier. I wouldn't think so, but I don't know your problems.

Shugojin
Sep 6, 2007

THE TAIL THAT BURNS TWICE AS BRIGHT...


In the following wall of ugly-rear end code, can someone PLEASE tell me what the s->f in the heavily indented lines is referring to.

I know that it means "member f of object s" but I'm not sure what that means in the rest of it. I have also indented and put !!! !!! around the line where I think the object s is defined but I'm lost as to what f is if that is the case.

And the entirety of this also belongs in the "Code that makes you want to laugh and/or cry" thread :suicide:

code:
#include </home/hank/include/gsl/gsl_vector.h>
#include </home/hank/include/gsl/gsl_blas.h>
#include </home/hank/include/gsl/gsl_multifit_nlin.h>

// this is a hack to make a numerical recipies compartible mrqmin function using gsl

struct data {
  int n;
  int np;
  float *p;
  float *dy;
  float *x;
  float *y;
  float *s;
  void (*f)(float, float *, float *, float *, int);
};

int fitfun_f(const gsl_vector *x, void *data, gsl_vector *f) {
  struct data *d = (struct data *) data;
  int n = d->n;
  int i;
  float y;
  for (i = 0; i < d->np; i++) {
    d->p[i] = (float) gsl_vector_get(x, i);
  }
  for (i = 0; i < n; i++) {
    //pointer d->p decrement 1, to accomodate the starting index of 1 in the gauss and lorentz functions.
    d->f(d->x[i], d->p-1, &y, NULL, d->np);
    gsl_vector_set(f, i, (y-d->y[i])/d->s[i]);

  }
  return GSL_SUCCESS;
}

int fitfun_df(const gsl_vector *x, void *data, gsl_matrix *J) {
  struct data *d = (struct data *) data;
   
  int n = d->n;
  int i, j;
  float y;
  for (i = 0; i < d->np; i++) {
    d->p[i] = (float) gsl_vector_get(x, i);
  }
  for (i = 0; i < n; i++) {
    //pointers x->data and d->dy decrement 1, to accomodate the starting index of 1 in the gauss and lorentz functions.
    d->f(d->x[i], d->p-1, &y, d->dy-1, d->np);
    for (j = 0; j < d->np; j++) {
      gsl_matrix_set(J, i, j, d->dy[j]/d->s[i]);
    }
  }
  return GSL_SUCCESS;
}

int fitfun_fdf(const gsl_vector *x, void *data, gsl_vector *f, gsl_matrix *J) {
  struct data *d = (struct data *) data;
  int n = d->n;
  int i, j;
  float y;
  for (i = 0; i < d->np; i++) {
    d->p[i] = (float) gsl_vector_get(x, i);
  }
  for (i = 0; i < n; i++) {
    //pointers x->data and d->dy decrement 1, to accomodate the starting index of 1 in the gauss and lorentz functions.
    d->f(d->x[i], d->p-1, &y, d->dy-1, d->np);
    gsl_vector_set(f, i, (y-d->y[i])/d->s[i]);
    for (j = 0; j < d->np; j++) {
      gsl_matrix_set(J, i, j, d->dy[j]/d->s[i]);
    }
  }
  return GSL_SUCCESS;
}    

const gsl_multifit_fdfsolver_type *T;
gsl_multifit_fdfsolver *s;
gsl_matrix *covar;
gsl_multifit_function_fdf f;
gsl_vector *p;
struct data d;

//does not make use of ia and alpha parameters.
void mrqmin(float *x,float *y,float *sig,int n,float *a,int *ia,int ma,
	    float **acovar,float **alpha,float *chisq,
	    void (*func)(float,float *,float *,float *,int),float *alambda) {
  int status, i, j;
  float ochisq;
  if (*alambda < 0.0) {
    d.n = n;
    d.np = ma;
    d.p = malloc(sizeof(float)*ma);
    d.dy = malloc(sizeof(float)*ma);
    d.x = x+1;
    d.y = y+1;
    d.s = sig+1;
    d.f = func;


    f.f = fitfun_f;
    f.df = fitfun_df;
    f.fdf = fitfun_fdf;
    f.n = n;
    f.p = ma;
    f.params = &d;

    covar = gsl_matrix_alloc(ma, ma);
    p = gsl_vector_alloc(ma);
    for (i = 1; i <= ma; i++) {
      gsl_vector_set(p, i-1, a[i]);
    }
    T = gsl_multifit_fdfsolver_lmsder;
    							!!! s = gsl_multifit_fdfsolver_alloc(T, n, ma); !!!
    gsl_multifit_fdfsolver_set(s, &f, p);

    *alambda = 0.001;
  }
  
  if (*alambda == 0) {
    gsl_multifit_covar(s->J, 0.0, covar);
    for (i = 1; i <= ma; i++) {
      for (j = 1; j <= ma; j++) {
	acovar[i][j] = (float)gsl_matrix_get(covar, i-1, j-1);
      }
    }
    free(d.p);
    free(d.dy);
    gsl_vector_free(p);
    gsl_matrix_free(covar);
    gsl_multifit_fdfsolver_free(s);
    return;
  }

  								ochisq = gsl_pow_2(gsl_blas_dnrm2(s->f));
  status = gsl_multifit_fdfsolver_iterate(s);
  								*chisq = gsl_pow_2(gsl_blas_dnrm2(s->f));
  if (*chisq <= ochisq) {
    *alambda *= 0.999;
  } else {
    *alambda *= 1.001;
  }
  for (i = 1; i <= ma; i++) {
    a[i] = (float)gsl_vector_get(s->x, i-1);
  }
}

roomforthetuna
Mar 22, 2005

I don't need to know anything about virii! My CUSTOM PROGRAM keeps me protected! It's not like they'll try to come in through the Internet or something!

Shugojin posted:

In the following wall of ugly-rear end code, can someone PLEASE tell me what the s->f in the heavily indented lines is referring to.
s->f appears to be a function [pointer], the value is presumably set in the call to multifit_fdfsolver_set.

Hammerite
Mar 9, 2007

And you don't remember what I said here, either, but it was pompous and stupid.
Jade Ear Joe

rjmccall posted:

Side note: I have no idea what you're trying to model here, and the signatures of these restrict methods make no sense to me. Is this a six-dimensional space? How do you restrict a coordinate range by a range in a different dimension, and what does the fixed point in the third dimension have to do with anything? It's hard, and probably fundamentally misguided of me, to suggest abstractions when I don't understand the problem space, but I'll keep trying anyway.

I've not really done anything to make it clear what this is about, and my posts have been pretty opaque from that perspective, so maybe I should briefly explain what I want to do in case this is a case of the X-Y problem :) It doesn't help that my thinking on this has been a little muddled, as well as me not being familiar with C++.

I'm not trying to model a six-dimensional space, it's actually only a two-dimensional one. I'm working with a hex grid. Now for a two-dimensional space you obviously only need two coordinates, but with a hex grid it is convenient for some purposes to introduce a third, linearly dependent coordinate. The neat way to do this is to have x + y + z = 0.



However, I also want to specify regions of the plane like this one:



This region cannot be described as the set of hexes satisfying x1 <= x <= x2, y1 <= y <= y2, z1 <= z <= z2 for some fixed integers x1, x2, y1, y2, z1, z2. However, write a = 2*x + y. Now we can describe this region by adding the requirement that a1 <= a <= a2 for some a1, a2. For symmetry, we add also two more coordinates b and c (which are honestly of more limited use given that the motivation for a was to describe oblong regions, but symmetry is nice) and we have a + b + c = 0. Call x, y and z the major coordinates and a, b, and c the minor coordinates. One issue is that pairing a major coordinate with a minor coordinate does not necessarily describe a hex; there is no hex such that a = 2, y = 1 (solving this we get x = 1/2).

I want to describe programmatically

- points (i.e. hexes) in the grid; these can be described by any two major coordinates, and we may as well systematically use x and y.
- lines in the grid; these can be described by one major coordinate
- "diagonal" lines in the grid - lines parallel to the a, b or c axes; these can be described by one minor coordinate


(perhaps an abuse of the words "orthogonal" and "diagonal")

- coordinate ranges; ranges of coordinate values like x1 <= x <= x2
- line segments; coordinate ranges through a point
- regions expressed as the intersections of 6 coordinate ranges, one in each of the 6 coordinates. It was while trying to write a function to calculate the area (number of hexes) of such a region that I started to question the way I was implementing this.

Perhaps I need to have MajorCoord and MinorCoord both inherit from Coord, have x, y, z inherit from the first and a, b, c inherit from the second, and then... something something... These ideas would be easier to juggle if I were practised in thinking in terms of C++'s type system. But I've always needed something I actually want to do to motivate me to learn a programming language, so there you go.

quote:

Anyway, you would obviously want these classes to share an implementation, instead of just copying them N times; among other things, that makes it a lot easier to abstract over them. I'd suggest:

I see this.

roomforthetuna
Mar 22, 2005

I don't need to know anything about virii! My CUSTOM PROGRAM keeps me protected! It's not like they'll try to come in through the Internet or something!
I always like to endorse doing hex grids with only two coordinates. If you rotate your grid to a particular angle, it can be represented as a square grid with six adjacencies to each square, north, northeast, east, south, southwest, west (ie. omitting the northwest-southeast diagonal). This makes it fit conveniently in computerese square arrays, albeit with maybe some wasted space in the corners depending on the shape of your whole grid.

You can convert such coordinates to screen/spatial coordinates using (screenx = hexsize*(x+0.5y)) [may be -0.5y I forget] and (screeny = hexsize*y)

Your square region of hexes would I think be represented by lowy<y<highy and lowa<a<higha where a=x+0.5y [or minus]

Or something like that, I'm not checking any of this for accuracy, but the gist works.

Lines can only be described by a single coordinate (in your system) if you're also assuming everything else is zero at one point on the line, and then you also need to know which of the six directions of line it is, so it's actually described by two things, not just by a single coordinate - having a "NorthEastLine" with X of 4 doesn't make the line described by just an X of 4.

I'm inclined to think you're trying to overengineer that aspect, and would be better off just representing your lines as three or four values, a two dimensional point on the line and a direction. (Three or four because the direction could be an index into a list of the few possible directions, or it could be a [1,0], depending which suits your goals better.) You could also represent it as a direction index and a 'perpendicular' direction distance, where perpendicular is defined as "one clockwise", but that seems likely to be more confusing and more difficult to work with than just describing lines more explicitly.

'Diagonal' directions would effectively be like [2,1] since that would be the step to the hexes they actually pass through. They jump in twos because they run the line between the adjacent hexes.

quote:

with a hex grid it is convenient for some purposes to introduce a third, linearly dependent coordinate
I have never seen it be more convenient for any purpose, but I would like to hear one! I've seen it done this way and it's always made an unnecessary, complicated, inefficient mess compared to the six-adjacency square grid method.

That Turkey Story
Mar 30, 2003

roomforthetuna posted:

I always like to endorse doing hex grids with only two coordinates. If you rotate your grid to a particular angle, it can be represented as a square grid with six adjacencies to each square, north, northeast, east, south, southwest, west (ie. omitting the northwest-southeast diagonal). This makes it fit conveniently in computerese square arrays, albeit with maybe some wasted space in the corners depending on the shape of your whole grid.
Seconding this. I'm not sure I see what your gains are from implementing a complicated coordinate system. Also, I'm not sure exactly what you're using this for, but before you start implementing algorithms for your hex tiles, look into the Boost Graph Library. It shouldn't be too difficult to make an implicit hex graph similar to Boost's grid_graph. Once you do that, you'll get things like path finding pretty much for free. I'm sort of surprised they don't have an implicit hex graph and triangle graph already.

Shugojin
Sep 6, 2007

THE TAIL THAT BURNS TWICE AS BRIGHT...


I wanna say thanks for your help guys I finally got that poo poo done. Ugh.

(gently caress Numerical Recipes in C right in its ported from Fortran rear end)

Deus Rex
Mar 5, 2005

I'm writing a raytracer for my first project in C++. For my primitive shapes, I've done something like this:

C++ code:
typedef vmath::Vector3d Point;
typedef vmath::Vector3d Vector;

// ...

class Plane
{
private:
    Vector normal;
    Point point;
public:
    Plane(Vector normal, Point point) : normal(normal), point(point) { } ;
}
(vmath library is here)

but since typedef only defines a type alias, a programmer could easily make an error like this:

C++ code:
Plane foo(Point(0,0,0), Vector(0,1,0));
I'd rather be able to catch this sort of error (semantic mixup of Points and Vectors) at compile time. Is there a standard way to do this? Maybe create empty Point and Vector classes which inherit from vmath::Vector3d?

Deus Rex fucked around with this message at 03:36 on Mar 9, 2013

Sarcophallus
Jun 12, 2011

by Lowtax

Deus Rex posted:

I'm writing a raytracer for my first project in C++. For my primitive shapes, I've done something like this:

C++ code:
typedef vmath::Vector3d Point;
typedef vmath::Vector3d Vector;

// ...

class Plane
{
private:
    Vector normal;
    Point point;
public:
    Plane(Vector normal, Point point) : normal(normal), point(point) { } ;
}
(vmath library is here)

but since typedef only defines a type alias, a programmer could easily make an error like this:

C++ code:
Plane foo(Point(0,0,0), Vector(0,1,0));
I'd rather be able to catch this sort of error (semantic mixup of Points and Vectors) at compile time. Is there a standard way to do this? Maybe create empty Point and Vector classes which inherit from vmath::Vector3d?

A quick google search led me to this, but I will readily admit that I've never used those and have no idea how you'd use them for this case.

An alternative is to pick one of your typedefs, and instead of doing the typedef for it, put it in a wrapper class. It's a bit round-a-bout, obviously, but it would at least generate an error for you if you mix types.

That Turkey Story
Mar 30, 2003

Deus Rex posted:

I'm writing a raytracer for my first project in C++. For my primitive shapes, I've done something like this:

Before you go any further, a plane can be defined by a normal and a scalar, but you are storing it as a normal and a point on the plane. It might be convenient to initialize a plane by providing a point and a normal, but it's not going to be very useful after that. Just calculate D and store that with the normal.

Deus Rex posted:

I'd rather be able to catch this sort of error (semantic mixup of Points and Vectors) at compile time. Is there a standard way to do this? Maybe create empty Point and Vector classes which inherit from vmath::Vector3d?
There are more differences as well, which makes a typedef insufficient anyway if you want to avoid misuse. Direct adding of points, for instance, doesn't make sense, though addition of a point and a vector does, and subtraction of two points does, etc. It's up to you how you want to go about representing things, obviously, but in practice, most people just use all vectors and refer to points as "position vectors." I don't like that, but in practice it's simpler and not too many libraries bother separating points from vectors. If you decide to go through and make a point type, either store a vector as a datamember or use private inheritance. Do not use public inheritance.

Adbot
ADBOT LOVES YOU

Jabor
Jul 16, 2010

#1 Loser at SpaceChem

That Turkey Story posted:

If you decide to go through and make a point type, either store a vector as a datamember or use private inheritance. Do not use public inheritance.

Public inheritance could be useful for interacting with libraries that just represent a Point as a vec3 or whatever, while still getting the benefits of type-safety in your own code. Assuming you can create all your Points yourself, that is - otherwise you'd have to wrap any library functions that return a Point, and then you might as well wrap everything in type-safe versions.

  • 1
  • 2
  • 3
  • 4
  • 5
  • Post
  • Reply