ISO_C_BINDING and arrays, optional arguments, etc.
vincenzo
working on my GTK+ / Fortran 2003 binding, I found among the thousand
functions (around 4600) of GTK+ one or two hundreds C prototypes of
the following kinds:
void gtk_signal_emit (GtkObject *object,guint signal_id,...);
void gtk_curve_get_vector (GtkCurve *curve,int veclen, gfloat
vector[]);
GtkWidget* gtk_clist_new_with_titles (gint columns,gchar *titles[]);
void gdk_window_invalidate_maybe_recurse (GdkWindow *window, const
GdkRegion *region, gboolean (*child_func) (GdkWindow *, gpointer),
gpointer user_data);
Have you some ideas to write interfaces for those functions using
ISO_C_BINDING ?
Vincent
James Van Buskirk
news:e6756b78-6e2c-411c...@r19g2000prm.googlegroups.com...
I suppose you're talking about GtkObject *Object, right? I'm still
not sure what you mean. One possibility is that GtkOject is an
opaque type where the user only sees pointers to it. In that case
you can either convert dummy arguments to TYPE(C_PTR),VALUE :: Object
and function results from GtkOject *fun() to TYPE(C_PTR) fun(), or
create an alias for C_PTR:
module alias_mod
use ISO_C_BINDING,ONLY: C_PTR ! Bug in older ifort requires C_PTR=>C_PTR
use ISO_C_BINDING,ONLY: GtkObject=>C_PTR
end module alias_mod
module main_mod
use alias_mod ! Now USErs of main_mod see GtkObject as an alias for C_PTR
end module main_mod
And now you can use GtkObject as you would have used C_PTR previously.
You still have the problem that a dummy argument could be declared
as GtkObject **pObject, in which case you have to check the
documentation to see whether the function expects an array of
GtkObject, a pointer to a pointer so that it can perhaps allocate
memory for a GtkObject and return a pointer to that memory through
pObject, or if he really wants a TYPE(C_PTR) by value that
dereferences to a pointer to GtkObject, or if he is just asking
for a reference to a pointer to GtkObject when passing the pointer
by value would have sufficed because it's INTENT(IN).
If you actually have the specification for GtkObject, then
obviously you want to write it out in Fortran. When you get
anything passed as type *t, if it's const type *t, then it
probably is an INTENT(IN) array, although I have seen a case in
wglext.h where it's an INTENT(IN) scalar -- you have to read the
documentation to make sure. If it's not const, than you have a
big problem because it may be an array of type, a scalar of type
passed by reference, or a TYPE(C_PTR) passed by value that is
required, and you have to read the documentation to make the
determination. That's where I am in glext.h: it took me about a day
to write a Fortran program that could parse all 1933 function
declarations and turn them into good Fortran like they should have
been written in in the first place, and the about 1000 instances of
pointers to constants all seemed to require INTENT(IN) arrays or raw
C_PTR by value, but the 400 or so non-const pointer dummies have to
be checked against the documentation on a case-by-case basis. I
got about 1/3 of the way through the task a week ago but I've been
fussing around with the GridView example from Norman Lawrence's
book since then. I'm hoping that when I get back to glext.h it
will only take about a day to finish, but the documentation is so
sparse that I don't know how anyone is ever supposed to try and
invoke some of those functions.
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
James Van Buskirk
news:iifb4u$rem$1...@news.eternal-september.org...
> "vincenzo" <vincent...@libertysurf.fr> wrote in message
> news:e6756b78-6e2c-411c...@r19g2000prm.googlegroups.com...
>> void gdk_window_invalidate_maybe_recurse (GdkWindow *window, const
>> GdkRegion *region, gboolean (*child_func) (GdkWindow *, gpointer),
>> gpointer user_data);
Oh, I didn't notice the procedure argument there. Write out a
prototype:
public PTRCHILD_FUNCPROC
abstract interface
function PTRCHILD_FUNCPROC(gpointer_) bind(C)
import
implicit none
!GCC$ ATTRIBUTES STDCALL :: PTRCHILD_FUNCPROC
type(GdkWindow) gpointer_ ! Or type(C_PTR),value :: gpointer_ if
opaque
integer(gboolean) PTRCHILD_FUNCPROC
end function PTRCHILD_FUNCPROC
end interface
Then you can write out an interface for the actual subroutine:
public gdk_window_invalidate_maybe_recurse
interface
subroutine gdk_window_invalidate_maybe_recurse(window,region, &
child_func,user_data) &
bind(C,name='gdk_window_invalidate_maybe_recurse')
import
implicit none
!GCC$ ATTRIBUTES STDCALL :: gdk_window_invalidate_maybe_recurse
type(GdkWindow) window
type(GdkRegion),intent(in) :: region
procedure(PTRCHILD_FUNCPROC) child_func
type(gpointer),value :: user_data
end subroutine gdk_window_invalidate_maybe_recurse
end interface
James Van Buskirk
> !GCC$ ATTRIBUTES STDCALL :: gdk_window_invalidate_maybe_recurse
I assumed that, as is usually the case, GTK+ uses STDCALL as the
calling convention for 32-bit Windows. I actually tried your
mandelbrot example today and found that this is not the case, so
the above Directive-Enhanced Compilation statement would be
incorrect. Also I observed the it took 25 seconds to render the
image in 64-bit Windows, but 3 minutes in 32-bit Windows!
James Van Buskirk
news:e6756b78-6e2c-411c...@r19g2000prm.googlegroups.com...
> Have you some ideas to write interfaces for those functions using
> ISO_C_BINDING ?
Oh, another thing I noticed was that some functions were using
LOGICAL(C_BOOL). Never do this, because no Fortran variable type
is interoperable with _Bool. It's a mistake for gfortran not to
assign C_BOOL = -1; using C_BOOL is just a recipe for disaster.
Start with gboolean = gint = int has storage size 4 bytes, but
_Bool has storage size 1 byte. I am also surprised to see no
C_INTPTR_T or C_SIZE_T declarations. Did none of the functions
contains gsize or gintptr types?
FX
> _Bool.
Would you care to back this affirmation up? Simple tests show that it
works nicely:
$ cat a.c
#include <stdio.h>
void foo (_Bool b)
{
printf ("Boolean is %s\n", b ? "true" : "false");
}
$ cat a.f90
use iso_c_binding
interface
subroutine foo(b) bind(c)
use iso_c_binding
logical(c_bool), value :: b
end subroutine
end interface
call foo(.true._c_bool)
call foo(.false._c_bool)
end
$ gfortran a.f90 a.c && ./a.out
Boolean is true
Boolean is false
> Start with gboolean = gint = int has storage size 4 bytes, but _Bool
> has storage size 1 byte.
So, what you really mean is that gboolean is not _Bool. That's quite
different from your overbroad assertion above.
glen herrmannsfeldt
(snip)
> Oh, another thing I noticed was that some functions were using
> LOGICAL(C_BOOL). Never do this, because no Fortran variable type
> is interoperable with _Bool. It's a mistake for gfortran not to
> assign C_BOOL = -1; using C_BOOL is just a recipe for disaster.
> Start with gboolean = gint = int has storage size 4 bytes, but
> _Bool has storage size 1 byte. I am also surprised to see no
> C_INTPTR_T or C_SIZE_T declarations. Did none of the functions
> contains gsize or gintptr types?
Who says that _Bool is one byte?
As I understand C interoperability, there is a corresponding
compiler which has the appropriate types and sizes.
There seem to be web pages indicating that gcc and gfortran do
have the appropriate sizes and values to make it work.
gcc claims that it does have a LOGICAL(KIND=C_BOOL) that
is interoperable with _Bool. Also that LOGICAL(KIND=N) variables
have the same size is INTEGER(KIND=N) variables.
Of course it won't work on a compiler that has C_BOOL .EQ. -1.
-- glen
vincenzo
> Also I observed the it took 25 seconds to render the
> image in 64-bit Windows, but 3 minutes in 32-bit Windows!
I had the same problem on an Intel Atom processor and Windows 7 +
gfortran 4.6. I removed all the GTK instructions to have just pure
fortran computation and the problem was the same. I tried also on
Windows XP + gfortran on a Pentium 4: same problem. But no problem
with Compaq Visual Fortran 6.5 under Windows.
So I think the problem comes either from my program or from gfortran
4.6 but not from GTK.
Vincent
vincenzo
> I suppose you're talking about GtkObject *Object, right? I'm still
> not sure what you mean.
void gtk_signal_emit (GtkObject *object,guint signal_id,...);
In fact, no, I was talking about ",..." which I think means that there
is optional arguments. Is it possible to deal with that in a Fortran
interface ?
But I am interested by your work on glext.h because I am also working
on a parsing program (in Python) to automatically write the
interfaces. With a little more work, I think it is possible to
automatically interface 90% of the GTK+ functions.
Concerning GtkWidget* gtk_clist_new_with_titles (gint columns,gchar
*titles[]);
my problem is with the array of C strings (which I usually interfaced
as an array of characters).
Vincent
Reinhold Bader
> On 3 fév, 23:44, "James Van Buskirk" <not_va...@comcast.net> wrote:
>> I suppose you're talking about GtkObject *Object, right? I'm still
>> not sure what you mean.
>
> void gtk_signal_emit (GtkObject *object,guint signal_id,...);
> In fact, no, I was talking about ",..." which I think means that there
> is optional arguments. Is it possible to deal with that in a Fortran
> interface ?
There are no provisions in the Fortran standard for directly interoperating with
a C vararg interface. However, you can write C wrappers to deal with
specific cases of signature combinations, and possibly a Fortran
generic which sells the various needed combinations under a single name.
Regards
Reinhold
nm...@cam.ac.uk
FX <cou...@alussinan.org> wrote:
>> Never do this, because no Fortran variable type is interoperable with
>> _Bool.
>
>Would you care to back this affirmation up? Simple tests show that it
>works nicely:
In trivial cases only. If you read the standards more carefully,
you will discover that they are semantically incompatible. Yes,
I know what can be done safely, but I would hate to teach it to
anyone without extensive Fortran, C and interoperability (and
preferably implementation) experience.
A lot of the Fortran interoperability section is a serious error,
because it does not have well-defined properties.
Regards,
Nick Maclaren.
James Van Buskirk
> In article <iiggi0$5qf$1...@news.eternal-september.org>,
> FX <cou...@alussinan.org> wrote:
>>> Never do this, because no Fortran variable type is interoperable with
>>> _Bool.
>>Would you care to back this affirmation up? Simple tests show that it
>>works nicely:
> In trivial cases only. If you read the standards more carefully,
> you will discover that they are semantically incompatible. Yes,
> I know what can be done safely, but I would hate to teach it to
> anyone without extensive Fortran, C and interoperability (and
> preferably implementation) experience.
Last couple of threads where I drew up some examples I just got
yelled at. In the case of gcc, the C compiler considers zero to
be .FALSE. and everything else to be .TRUE., but the Fortran
compiler considers zero to be .FALSE., one to be .TRUE., and
intentionally abuses you if you put any other value into a LOGICAL
variable. Thus code that read an _Bool value and tells the
Fortran compiler that it's a LOGICAL(C_BOOL) variable may encounter
any kind of mayhem when the optimizer starts composing syllogisms.
On the other hand, reading an _Bool value as an INTEGER(C_INT8_T)
and then testing whether it's nonzero will be safe in gcc. The
way to fix this is for the C compiler to to consider the LSB of
an _Bool variable to be the only significant bit, as the standard
allows, and then to follow suit with the Fortran compiler, but it's
so ingrained in the mindset of the C programmer that zero is .FALSE.
and everything else is .TRUE. that this situation will never happen
no matter how error-prone it is. Time for the yelling, I suppose.
But if you do that I will just have to put another Yes album on
and turn the volume way up so that I can't hear you any more.
> A lot of the Fortran interoperability section is a serious error,
> because it does not have well-defined properties.
Could you be kind enough to point out maybe a couple of the most
glaring cases?
James Van Buskirk
> On 3 fév, 23:44, "James Van Buskirk" <not_va...@comcast.net> wrote:
> > I suppose you're talking about GtkObject *Object, right? I'm still
> > not sure what you mean.
> void gtk_signal_emit (GtkObject *object,guint signal_id,...);
> In fact, no, I was talking about ",..." which I think means that there
> is optional arguments. Is it possible to deal with that in a Fortran
> interface ?
Nope. In this case you're dead meat. This syntax is specifically
called out in the standard as uninteroperable. What would ifort do?
"ATTRIBUTES VARYING
The ATTRIBUTES directive option VARYING allows a variable number of
calling arguments. It takes the following form:
cDEC$ ATTRIBUTES VARYING :: var [, var] ...
c
Is one of the following: C (or c), !, or *. (See Syntax Rules for
Compiler Directives.)
var
Is the name of a variable.
Either the first argument must be a number indicating how many
arguments to process, or the last argument must be a special marker
(such as -1) indicating it is the final argument. The sequence of
the arguments, and types and kinds must be compatible with the called
procedure.
If VARYING is specified, the C option must also be specified."
So it doesn't sound promising even on a compiler that implements
varying arguments. The thing to do if someone actually decides to
torture you with this stuff is to write out an interface block for
the specific instance you are being made to suffer with. Adn hope
there is no quirk about varying arguments that makes it fail on any
platform you are interested in.
> But I am interested by your work on glext.h because I am also working
> on a parsing program (in Python) to automatically write the
> interfaces. With a little more work, I think it is possible to
> automatically interface 90% of the GTK+ functions.
If it took a sophisticated parser to do the job, it wouldn't get done,
at least not by me. I just started out with Fortran code that would
handle a couple of simple cases and harshly reject anything that
didn't fit within my limited expectations. When the first exception
occurred, I wrote to code handle it and the ran the parser again.
Repeating this process for several hours finally resulted in a program
that could fight its way it through the file. The process was helped
greatly by code written in a uniform style (I mean the C code, not
mine!)
> Concerning GtkWidget* gtk_clist_new_with_titles (gint columns,gchar
> *titles[]);
> my problem is with the array of C strings (which I usually interfaced
> as an array of characters).
Now that one is easy:
type(C_PTR) titles(*)
nm...@cam.ac.uk
>> because it does not have well-defined properties.
>
>Could you be kind enough to point out maybe a couple of the most
>glaring cases?
Well, this one - which also includes C_SIZE_T.
C_ASSOCIATED. What the HELL does associated mean in a C context?
C has null pointers, addresses of lvalues, valid addresses but
with no lvalue, valid pointer values that are not addresses, and
undefined values. Which? And what does C_ASSOCIATED(x,y) do if
x is a subobject of y? And what IS an entity in this context,
anyway?
C_F_POINTER. What are the rules for type compatibility, both when
the values will be inspected, and when only the storage is used?
C and Fortran have WILDLY different rules, and the former is a
spectacular mess. For example, what on earth does the 'type of
the entity' mean for something pointed to by void * (as C_PTR
indicates)?
Regards,
Nick Maclaren.
JB
[*sigh*, not this canard again]
> In the case of gcc, the C compiler considers zero to
> be .FALSE. and everything else to be .TRUE.
No, this is not correct for _Bool. You might want to ponder the
assembler output from
_Bool neg(_Bool b)
{
return !b;
}
and compare it to the asm generated for
int neg(int b)
{
return !b;
}
and for kicks, you might even check the equivalent
logical function neg(b)
logical :: b
neg = .not. b
end function neg
with and without kind=C_BOOL.
>, but the Fortran
> compiler considers zero to be .FALSE., one to be .TRUE., and
> intentionally abuses you if you put any other value into a LOGICAL
> variable.
Same for _Bool. Yes, really.
So yes, this does mean that logical operations on _Bool variables can
be implemented differently than the same operations on other primitive
types. Perhaps one reason why C programmers haven't been particularly
enthusiastic about _Bool.. :)
> The
> way to fix this is for the C compiler to to consider the LSB of
> an _Bool variable to be the only significant bit, as the standard
> allows,
The standard might allow that, yes, but it sounds like a somewhat
impractical choice..
--
JB
James Van Buskirk
news:slrnikopa...@kosh.hut.fi...
[snip]
I don't have time to get you up to speed from zero. Go back to the
last couple of threads about this and refute specific quotations I
made there from both standard documents and also refute the
examples I posted in both C and Fortran if that's what you want to
do, but I just don't have time to reestablish the context of those
previous discussions. I'm falling behind on projects that are
enjoyable right now.
James Van Buskirk
I just had a thought that the difference may be SSE2, so I tried:
gfortran -Wall -Wno-unused-dummy-argument -O2 -march=core2 -mfpmath=sse
mandelbrot.f90 -omandelbrot gtk.o -LC:\GTK32\lib -lgtk-win32-2.0.dll
-lgobject-2.0.dll -lgdk-win32-2.0.dll
But it didn't help. Also that -Wall turned up some tabs and unused
arguments that you might want to clean up. If you want to see fast
mandelbrot stuff, check out:
http://board.flatassembler.net/topic.php?t=5122
JB
> last couple of threads about this and refute specific quotations I
> made there from both standard documents and also refute the
> examples I posted in both C and Fortran if that's what you want to
> do,
Hmm, why should I do that? I'm not the one wanting to change the
status quo, and as I recall, I didn't find your arguments back then
particularly convincing. So no, I'm not going to spend time on a wild
goose chase trying to prove you wrong in a way that satisfies your
highness.
In any case, if you want to change how gcc handles _Bool variables, my
opinion doesn't particularly matter; rather, you need to convince the
gcc C frontend maintainers that their interpretation of the C standard
is incorrect. If you manage that, I find it likely that gfortran would
follow, as the frontends share the underlying machinery for handling
boolean types, and the Fortran standard allows the internal
representation of logical variables to be whatever the compiler
pleases.
--
JB
nm...@cam.ac.uk
>On 2011-02-04, James Van Buskirk <not_...@comcast.net> wrote:
>> Go back to the
>> last couple of threads about this and refute specific quotations I
>> made there from both standard documents and also refute the
>> examples I posted in both C and Fortran if that's what you want to
>> do,
>
>opinion doesn't particularly matter; rather, you need to convince the
>gcc C frontend maintainers that their interpretation of the C standard
>is incorrect. If you manage that, I find it likely that gfortran would
>follow, as the frontends share the underlying machinery for handling
>boolean types, and the Fortran standard allows the internal
>representation of logical variables to be whatever the compiler
>pleases.
If you are claiming that C99 _Bool has compatible semantics to
Fortran LOGICAL, either in the standard or gcc, you are talking
nonsense. I accept that you can't expose the differences without
explicit or implicit casting, but that is a fully-specified part
of the C standard.
There are also subtle differences between C99 _Bool and C++ bool,
and a huge proportion of people using the interoperability are
actually interfacing into the C subset of C++. That's a gotcha.
I support James Van Buskirk here - the situation is just too full
of subtle gotchas to be safe for the ordinary programmer. There
IS a path through the minefield, but you had better know all of
the Fortran, C and C++ standards very well to be sure of following
it correctly.
Regards,
Nick Maclaren.
JB
> In article <slrnikqar...@kosh.hut.fi>, JB <f...@bar.invalid> wrote:
>>In any case, if you want to change how gcc handles _Bool variables, my
>>opinion doesn't particularly matter; rather, you need to convince the
>>gcc C frontend maintainers that their interpretation of the C standard
>>is incorrect. If you manage that, I find it likely that gfortran would
>>follow, as the frontends share the underlying machinery for handling
>>boolean types, and the Fortran standard allows the internal
>>representation of logical variables to be whatever the compiler
>>pleases.
>
> If you are claiming that C99 _Bool has compatible semantics to
> Fortran LOGICAL, either in the standard or gcc, you are talking
> nonsense.
No, I'm not. My argument is basically that
1) Gfortran's implementation of LOGICAL conforms to the Fortran
standard. In particular, cute type-punning tricks are outside the
scope of the standard.
2) As it happens, the internal representation of LOGICAL variables
chosen by gfortran is the same as the one gcc uses for _Bool, allowing
both frontends to share some of the underlying machinery (e.g. code
generation for logical expressions). This does not exclude that the
frontends can implement different semantics for the variables.
3) Looking at the asm generated by gcc, it seems the C frontend
maintainers believe that poking at the internals of a _Bool variable
such that the corresponding integer value is neither 0 or 1, is not
allowed. My own reading of C99, which I do trust much less than that
of the C frontend maintainers, some of which AFAIK have contributed to
the standards effort, seems to confirm this.
> There are also subtle differences between C99 _Bool and C++ bool,
> and a huge proportion of people using the interoperability are
> actually interfacing into the C subset of C++. That's a gotcha.
I haven't looked at this, particularly the C++ standard, in sufficient
detail to be able to comment on this. Sorry.
> I support James Van Buskirk here - the situation is just too full
> of subtle gotchas to be safe for the ordinary programmer. There
> IS a path through the minefield, but you had better know all of
> the Fortran, C and C++ standards very well to be sure of following
> it correctly.
So what should a Fortran compiler do then, in your opinion? I suspect
setting C_BOOL=-1 and telling users "sorry, some corner case is too
complicated for you to handle so we won't support this at all" would
not go down too well.
--
JB
nm...@cam.ac.uk
>>
>> If you are claiming that C99 _Bool has compatible semantics to
>> Fortran LOGICAL, either in the standard or gcc, you are talking
>> nonsense.
>
>No, I'm not. My argument is basically that
>
>1) Gfortran's implementation of LOGICAL conforms to the Fortran
>standard. In particular, cute type-punning tricks are outside the
>scope of the standard.
Agreed.
>2) As it happens, the internal representation of LOGICAL variables
>chosen by gfortran is the same as the one gcc uses for _Bool, ...
Quite reasonably.
>3) Looking at the asm generated by gcc, it seems the C frontend
>maintainers believe that poking at the internals of a _Bool variable
>such that the corresponding integer value is neither 0 or 1, is not
>allowed. My own reading of C99, which I do trust much less than that
>of the C frontend maintainers, some of which AFAIK have contributed to
>the standards effort, seems to confirm this.
You are entirely right to distrust C99, but you are wrong about
what it says; there was a heated debate on this, which was resolved
in the usual way, by denying that there is a problem. The ambiguity
is to do with initialisation, but the following is clearly allowed:
_Bool x;
*(char *)x = 8;
God help us all, but the assignment "x = x;" will then set x to 1!
>So what should a Fortran compiler do then, in your opinion? I suspect
>setting C_BOOL=-1 and telling users "sorry, some corner case is too
>complicated for you to handle so we won't support this at all" would
>not go down too well.
It depends on whether it regards RAS or the opinions of the less
clued-up parts of the C community as more important. C99 has got
the thumbs-down from a large proportion of the C community, and
probably the majority of those that interface with Fortran.
Regards,
Nick Maclaren.
FX
Yet, it's very widely used.
--
FX
James Van Buskirk
> If you are claiming that C99 _Bool has compatible semantics to
> Fortran LOGICAL, either in the standard or gcc, you are talking
> nonsense. I accept that you can't expose the differences without
> explicit or implicit casting, but that is a fully-specified part
> of the C standard.
Now, if you could provide an example where the C companion processor
and the Fortran processor both agree that the LSB of an _Bool or a
LOGICAL(C_BOOL) were the only significant bit and all the other bits
were copied on assignment to another _Bool or LOGICAL(C_BOOL) but
yet the semantics proved somehow incompatible or inconsistent, you
would be doing me a great service because I would then be able to
throw my hands up and pronounce the whole situation impossible and
not get drawn into discussions where so many posters seem to get
their kicks from abusing me.
nm...@cam.ac.uk
>> Fortran LOGICAL, either in the standard or gcc, you are talking
>> nonsense. I accept that you can't expose the differences without
>> explicit or implicit casting, but that is a fully-specified part
>> of the C standard.
>
>Now, if you could provide an example where the C companion processor
>and the Fortran processor both agree that the LSB of an _Bool or a
>LOGICAL(C_BOOL) were the only significant bit and all the other bits
>were copied on assignment to another _Bool or LOGICAL(C_BOOL) but
>yet the semantics proved somehow incompatible or inconsistent, you
>would be doing me a great service because I would then be able to
>throw my hands up and pronounce the whole situation impossible and
>not get drawn into discussions where so many posters seem to get
>their kicks from abusing me.
Indeed, but that is not the situation. Regrettably, like so much
of C (and some other aspects of Fortran/C interoperability), it is
unclear to the point of ambiguity and even inconsistency.
I believe that they CAN BE interoperable. Some combinations of
compilers will do that, and will set C_BOOL to -1, and some won't
but will set C_BOOL inappropriately. Similar remarks apply about
C++ bool and C99 _Bool.
I think that, only with the shoddiest of compiler combinations,
will the last occur. HOWEVER, the interoperability I am referring
to is with the SUBSET of C that makes sense in both Fortran and C++
terms - i.e. where Boolean objects are never updated via an lvalue
of any other type.
What I can also say is that it will be beyond the ability of most
support staff to get a vendor to accept a bug report in a case
where there is a non-trivial failure of interoperability, even
with that assumption. gcc/gfortran/g++ is probably an exception.
The point is that there is so much murk that it needs someone with
at least my level of knowledge of the standards and obstinacy to
argue against every attempt to weasel out.
For example, the integer promotions (and casts) preserve value,
and it is unclear whether a Fortran LOGICAL necessarily compares
equal to C true. Some compilers have had that fail, despite a
claim of (pre-standard) interoperability.
Regards,
Nick Maclaren.
James Van Buskirk
If you want to see the kind of progress I am making getting OpenGL
to work in gfortran on Windows, I have a new example at:
http://home.comcast.net/~kmbtib/Fortran_stuff/gridview.zip
Which has source code for the program and OpenGL32 interfaces as
well as pre-built 32- and 64-bit Windows executables. It's an
upgrade of the GridView example from Norman Lawrence's book:
http://www.elsevierdirect.com/product.jsp?isbn=1555582494
where the source code for the original and a 32-bit Windows
executable may be found as well. Let me know what you think.
baf
> "vincenzo"<vincent...@libertysurf.fr> wrote in message
> news:9e7b7e62-4eea-4450...@w36g2000vbi.googlegroups.com...
>
>
>>> Also I observed the it took 25 seconds to render the
>>> image in 64-bit Windows, but 3 minutes in 32-bit Windows!
>
>> I had the same problem on an Intel Atom processor and Windows 7 +
>> gfortran 4.6. I removed all the GTK instructions to have just pure
>> fortran computation and the problem was the same. I tried also on
>> Windows XP + gfortran on a Pentium 4: same problem. But no problem
>> with Compaq Visual Fortran 6.5 under Windows.
>
>> So I think the problem comes either from my program or from gfortran
>> 4.6 but not from GTK.
>
> If you want to see the kind of progress I am making getting OpenGL
> to work in gfortran on Windows, I have a new example at:
>
> http://home.comcast.net/~kmbtib/Fortran_stuff/gridview.zip
>
> Which has source code for the program and OpenGL32 interfaces as
> well as pre-built 32- and 64-bit Windows executables. It's an
> upgrade of the GridView example from Norman Lawrence's book:
>
> http://www.elsevierdirect.com/product.jsp?isbn=1555582494
>
> where the source code for the original and a 32-bit Windows
> executable may be found as well. Let me know what you think.
>
Wow! Really nice looking plot, but there sure is a lot of code required
to produce the figure!
James Van Buskirk
Since I got the gridview stuff done, I thought I would take a look
at the disassembly of your code. There is an invocation of _cabs
that may still invoke an x87 function even when -march=core2
-mfpmath=sse is in force. I don't know for sure, but I don't think
gfortran supplies more than one math library per installation.
Even if that is pure superstition on my part, making the following
change:
! do while ((k <= itermax) .and. (abs(z)<2d0))
do while ((k <= itermax) .and. (real(conjg(z)*z)<4d0))
took execution time from 3 minutes down to 7 seconds for a 32-bit
build.
vincenzo
> If you want to see the kind of progress I am making getting OpenGL
> to work in gfortran on Windows, I have a new example at:
>
> http://home.comcast.net/~kmbtib/Fortran_stuff/gridview.zip
>
> Which has source code for the program and OpenGL32 interfaces as
> well as pre-built 32- and 64-bit Windows executables. It's an
> upgrade of the GridView example from Norman Lawrence's book:
>
> http://www.elsevierdirect.com/product.jsp?isbn=1555582494
>
> where the source code for the original and a 32-bit Windows
> executable may be found as well. Let me know what you think.
Hello James, this is impressive ! (and the .exe file is even working
fine under Linux with Wine)
When I have time I will look at your interface Opengl32.f90. It can
inspire me for my gtk-fortran binding.In the devel branch I have now
around 8000 GTK+ functions (probably with some type errors ?). You can
see the advances of the project here:
https://github.com/jerryd/gtk-fortran/wiki
I saw in the GTK+ documentation that there are OpenGL GtkWidgets:
http://library.gnome.org/devel/gtkglext/stable/gtkglext-gtkglwidget.html
http://library.gnome.org/devel/gtkglext/stable/gtkglext-gdkglwindow.html
http://projects.gnome.org/gtkglext/
James Van Buskirk
> On 11 fév, 04:27, "James Van Buskirk" <not_va...@comcast.net> wrote:
> > If you want to see the kind of progress I am making getting OpenGL
> > to work in gfortran on Windows, I have a new example at:
> > http://home.comcast.net/~kmbtib/Fortran_stuff/gridview.zip
> > Which has source code for the program and OpenGL32 interfaces as
> > well as pre-built 32- and 64-bit Windows executables. It's an
> > upgrade of the GridView example from Norman Lawrence's book:
> > http://www.elsevierdirect.com/product.jsp?isbn=1555582494
> > where the source code for the original and a 32-bit Windows
> > executable may be found as well. Let me know what you think.
> Hello James, this is impressive ! (and the .exe file is even working
> fine under Linux with Wine)
I am curious about whether two of the more difficult features work
with Wine: printing (try File -> Print) and palette management (try
running gridview while the display adapter is set to 8 bit color
depth.)
> When I have time I will look at your interface Opengl32.f90. It can
> inspire me for my gtk-fortran binding.In the devel branch I have now
> around 8000 GTK+ functions (probably with some type errors ?). You can
> see the advances of the project here:
> https://github.com/jerryd/gtk-fortran/wiki
It took me a while to find where the new stuff is in your web page.
I'm trying to finish up the OpenGL extension stuff now.
> I saw in the GTK+ documentation that there are OpenGL GtkWidgets:
> http://library.gnome.org/devel/gtkglext/stable/gtkglext-gtkglwidget.html
> http://library.gnome.org/devel/gtkglext/stable/gtkglext-gdkglwindow.html
> http://projects.gnome.org/gtkglext/
Looking at the last link the project doesn't seem to have been
maintained in a couple of years and the page with the *.dll I
would need to make this work was 404 (in Italian!) so this looks
like another of those near misses. I do have a GLUT example on
my web page:
http://home.comcast.net/~kmbtib/Fortran_stuff/glut_test.zip
which should compile on any platform where you can find and
install freeglut along with gfortran, but GLUT lacks flexibility.
I don't know if it can render to the printer like you can with
Win32 API or if it can do palette management. Although I just
now tried running 12.exe from the glut_test.zip file above in
8-bit color depth and it looked quite deficient, but when I
switched to 32-bit color depth with the program still running
the colors came out looking completely wacky and beautiful.
You've just got to try it out to see what I mean. Running
another instance of 12.exe side by side with normal colors also
looks cool.
vincenzo
> Even if that is pure superstition on my part, making the following
> change:
>
> ! do while ((k <= itermax) .and. (abs(z)<2d0))
> do while ((k <= itermax) .and. (real(conjg(z)*z)<4d0))
>
> took execution time from 3 minutes down to 7 seconds for a 32-bit
> build.
Than you James. I confirm it is now faster under Windows 7 !
V.
vincenzo
> I am curious about whether two of the more difficult features work
> with Wine: printing (try File -> Print) and palette management (try
> running gridview while the display adapter is set to 8 bit color
> depth.)
Concerning Printing, it does nothing. I just obtain the following
error in console with wine 1.3.5:
fixme:commdlg:PrintDlgExA (0x44f140) not fully implemented
Concerning, 8 bit color depth, I did not succeed to write to xorg.conf
with my Nvidia settings and I prefer not insisting to keep my portrait
screen keep working...
But everything else seems to work with wine.
> It took me a while to find where the new stuff is in your web page.
OK, I will add some instructions on the Wiki.
Gordon Sande
The hint is in the change of the test value from 2.0d0 to 4.0d0. No
superstition involved.
ABS will use algorithms which go to some trouble to avoid overflows when
the value is greater than the square root of the maximum value.
Typically that involves finding the sup norm of the z (maximum absolute
value of the real or imaginary component) and then finding the
modificaion needed to adjust that to be to be the L_2 norm (a factor
between 1.0 and sqrt(2.0)).
The BLAS has examples of the work required. Look at the norms of complex
vectors.
Ron Shepard
Gordon Sande <Gordon...@gmail.com> wrote:
> On 2011-02-13 10:41:19 -0400, vincenzo said:
>
This is the so-called "Pythagorean Sum" problem which occurs in many
contexts including complex abs(). See the LAPACK routine DLAPY2 for
example,
http://www.netlib.org/lapack/explore-html/d5/d7c/dlapy2_8f_source.htm
l
which scales the factors in such a way to avoid overflow for large
values and then uses the SQRT() function. It can also be computed
iteratively without the SQRT(), but I think enough hardware supports
SQRT() these days so that the LAPACK approach is probably best for a
general routine. The Pythagorean sum is an expensive operation, and
math libraries, including LAPACK, often use extensively other scalar
and vector norms for complex arguments for which the overflow
problem is easier to solve and that do not require the sqrt(). For
example, instead of testing the above to make sure that Z is within
the circle of radius 2, the test might instead be changed to test
that it is within the square with sides 4.0 or 2.828 as appropriate,
or that ABS(REAL(Z))+ABS(IMAG(Z)) or MAX(ABS(REAL(Z)),ABS(IMAG(Z)))
is smaller than the appropriate value. These tests requires only
ABS() of the real and imaginary parts of Z, reduce or eliminate the
risk of overflow, and they are often faster even than the complex
multiplication used in the above replacement for the Pythagorean sum
test.
Avoiding overflow is very important on some machines. In VAX
hardware, for example, one of the 64-bit floating point
representations had the feature that X*X would overflow when X had a
magnitude of 1.0/EPS (or maybe just slightly larger than that
value). It was quite common when converging some value down to the
last couple of bits for some subexpression to overflow because of
that, so you constantly had to worry about scaling subexpressions
appropriately. In fact, I remember "fixing" one of these problems
once by changing the expression sqrt(x*x+y*y) to the nonstandard
cdabs(dcmplx(x,y)) to take advantage of the fact that the VAX
fortran intrinsic library was written correctly in this case. Of
course, I had to fix that "fix" the next time I tried to compile the
code on another machine.
So, complex ABS(Z) is an expensive operation, it can be avoided with
REAL(CONJG(Z)*Z) as above if you know that the multiplication will
never overflow, and if the convergence test can be generalized, it
can be replaced with even simpler and cheaper tests that do not risk
overflow in the first place.
$.02 -Ron Shepard
glen herrmannsfeldt
(big snip)
> So, complex ABS(Z) is an expensive operation, it can be avoided with
> REAL(CONJG(Z)*Z) as above if you know that the multiplication will
> never overflow, and if the convergence test can be generalized, it
> can be replaced with even simpler and cheaper tests that do not risk
> overflow in the first place.
Is that better than REAL(Z)**2+AIMAG(Z)**2?
It would seem that it would otherwise compute the useless cross term.
-- glen
James Van Buskirk
news:ij8sml$brg$1...@news.eternal-september.org...
> The hint is in the change of the test value from 2.0d0 to 4.0d0. No
> superstition involved.
> ABS will use algorithms which go to some trouble to avoid overflows when
> the value is greater than the square root of the maximum value.
> Typically that involves finding the sup norm of the z (maximum absolute
> value of the real or imaginary component) and then finding the
> modificaion needed to adjust that to be to be the L_2 norm (a factor
> between 1.0 and sqrt(2.0)).
> The BLAS has examples of the work required. Look at the norms of complex
> vectors.
OK, let's do so:
FUNCTION SCNRM2( N, X, INCX )
real(wp) scnrm2
* .. Scalar Arguments ..
INTEGER, intent(in) :: INCX, N
* .. Array Arguments ..
COMPLEX(wp), intent(in) :: X( * )
* ..
*
* SCNRM2 returns the euclidean norm of a vector via the function
* name, so that
*
* SCNRM2 := sqrt( conjg( x' )*x )
*
*
*
* -- This version written on 25-October-1982.
* Modified on 14-October-1993 to inline the call to CLASSQ.
* Sven Hammarling, Nag Ltd.
*
! modified 12/11/04, intents and real kind added - James Van Buskirk
*
* .. Parameters ..
REAL(wp) ONE , ZERO
PARAMETER ( ONE = 1.0E+0_wp, ZERO = 0.0E+0_wp )
* .. Local Scalars ..
INTEGER IX
REAL(wp) NORM, SCALE, SSQ, TEMP
* .. Intrinsic Functions ..
INTRINSIC ABS, AIMAG, REAL, SQRT
* ..
* .. Executable Statements ..
IF( N.LT.1 .OR. INCX.LT.1 )THEN
NORM = ZERO
ELSE
SCALE = ZERO
SSQ = ONE
* The following loop is equivalent to this call to the LAPACK
* auxiliary routine:
* CALL CLASSQ( N, X, INCX, SCALE, SSQ )
*
DO 10, IX = 1, 1 + ( N - 1 )*INCX, INCX
IF( REAL( X( IX ) ).NE.ZERO )THEN
TEMP = ABS( REAL( X( IX ) ) )
IF( SCALE.LT.TEMP )THEN
SSQ = ONE + SSQ*( SCALE/TEMP )**2
SCALE = TEMP
ELSE
SSQ = SSQ + ( TEMP/SCALE )**2
END IF
END IF
IF( AIMAG( X( IX ) ).NE.ZERO )THEN
TEMP = ABS( AIMAG( X( IX ) ) )
IF( SCALE.LT.TEMP )THEN
SSQ = ONE + SSQ*( SCALE/TEMP )**2
SCALE = TEMP
ELSE
SSQ = SSQ + ( TEMP/SCALE )**2
END IF
END IF
10 CONTINUE
NORM = SCALE * SQRT( SSQ )
END IF
*
SCNRM2 = NORM
RETURN
*
* End of SCNRM2.
*
END function scnrm2
Surprising that f08 doesn't allow NORM2 of a complex vector. I
suppose the compiler should know how to optimize NORM2(ABS(z)).
The above code seems rather wrongheaded. In my opinion, the loop
should process in 4 modes: initial mode, can't underflow mode,
underflow protect, and overflow protect.
In initial mode, we don't know whether we may have to apply a
scale factor to protect against underflow or overflow. In this
mode the absolute value of every real and imaginary part should
be checked. If one is so large as to create danger of overflow,
switch to overflow protect mode. Else if one is large enough that
values that would underflow on squaring switch to can't underflow
mode. Else if a value that would underflow on squaring is
encountered, switch to underflow protect mode.
In can't underflow mode, just check each part to see if we should
switch to overflow protect mode.
In underflow protect mode, check each part to see whether it is
so large that all previous work will not affect the result. If
so, discard al previous work and determine whether to switch to
can't underflow or overflow protect mode. If the value isn't
that large, multiply by a suitable power of 2 then accumulate
into the result.
In oveflow protect mode, multiply by a suitable power of 2 and
accumulate into the result.
At the end, we may have to multiply by a power of 2 if we were
in underflow or overflow protect mode. Lots of comparisons to
be sure, but you can't have that many mispredicted braches
over the course of a large scnrm2 calculation, and there are
no divisions. Also we never incur roundoff error due to the
correction factors because they are all powers of 2. Maybe
for radix-10 floating point or base 16 exponents we might have
to choose our scale factors a little differently to avoid
roundoff error due to scaling.
Here is a possible rewrite of cabs:
C:\gfortran\clf\csqrt>type cabs.f90
elemental function cabs_gen(z)
implicit none
complex(wp),intent(in),value :: z
real(wp),parameter :: x = 0
real(wp),parameter :: tolhi = nearest(sqrt(huge(x)/2),-1.0)
real(wp),parameter :: tollo = nearest(sqrt(tiny(x)/epsilon(x)),1.0)
real(wp),parameter :: mulhi = scale(1.0_wp,exponent(tolhi/huge(x))-1)
real(wp),parameter :: divhi = scale(1.0_wp,1-exponent(tolhi/huge(x)))
real(wp),parameter :: mullo = scale(1.0_wp,exponent(tollo/tiny(x))+1)
real(wp),parameter :: divlo = scale(1.0_wp,-1-exponent(tollo/tiny(x)))
real(wp) cabs_gen
if(abs(real(z)) > tolhi) then
cabs_gen = divhi*sqrt((mulhi*real(z))**2+(mulhi*aimag(z))**2)
else if(abs(real(z)) < tollo) then
if(abs(aimag(z)) < tollo) then
cabs_gen = divlo*sqrt((mullo*real(z))**2+(mullo*aimag(z))**2)
else if(abs(aimag(z)) > tolhi) then
cabs_gen = abs(aimag(z))
else
cabs_gen = sqrt(real(z)**2+aimag(z)**2)
end if
else
if(abs(aimag(z)) > tolhi) then
cabs_gen = divhi*sqrt((mulhi*real(z))**2+(mulhi*aimag(z))**2)
else
cabs_gen = sqrt(real(z)**2+aimag(z)**2)
end if
end if
end function cabs_gen
C:\gfortran\clf\csqrt>type cabs_test.f90
module m4
use ISO_FORTRAN_ENV
implicit none
private
integer, parameter :: wp = REAL_KINDS(1)
public cabs
interface cabs
module procedure cabs_gen
end interface cabs
contains
include 'cabs.f90'
end module m4
module m8
use ISO_FORTRAN_ENV
implicit none
private
integer, parameter :: wp = REAL_KINDS(2)
public cabs
interface cabs
module procedure cabs_gen
end interface cabs
contains
include 'cabs.f90'
end module m8
module m10
use ISO_FORTRAN_ENV
implicit none
private
integer, parameter :: wp = REAL_KINDS(3)
public cabs
interface cabs
module procedure cabs_gen
end interface cabs
contains
include 'cabs.f90'
end module m10
module m16
use ISO_FORTRAN_ENV
implicit none
private
integer, parameter :: wp = REAL_KINDS(4)
public cabs
interface cabs
module procedure cabs_gen
end interface cabs
contains
include 'cabs.f90'
end module m16
program testme
use ISO_FORTRAN_ENV
use m4
use m8
use m10
use m16
implicit none
complex(REAL_KINDS(1)) x4
complex(REAL_KINDS(2)) x8
complex(REAL_KINDS(3)) x10
complex(REAL_KINDS(4)) x16
x4 = cmplx(huge(real(x4))/10,huge(real(x4))/10,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(huge(real(x4))/10,1,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(huge(real(x4))/10,tiny(real(x4))*10,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(1,huge(real(x4))/10,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(tiny(real(x4))*10,huge(real(x4))/10,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(1,1,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(1,tiny(real(x4))*10,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(tiny(real(x4))*10,1,kind(x4))
write(*,*) x4, cabs(x4)
x4 = cmplx(tiny(real(x4))*10,tiny(real(x4))*10,kind(x4))
write(*,*) x4, cabs(x4)
write(*,*)
x8 = cmplx(huge(real(x8))/10,huge(real(x8))/10,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(huge(real(x8))/10,1,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(huge(real(x8))/10,tiny(real(x8))*10,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(1,huge(real(x8))/10,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(tiny(real(x8))*10,huge(real(x8))/10,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(1,1,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(1,tiny(real(x8))*10,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(tiny(real(x8))*10,1,kind(x8))
write(*,*) x8, cabs(x8)
x8 = cmplx(tiny(real(x8))*10,tiny(real(x8))*10,kind(x8))
write(*,*) x8, cabs(x8)
write(*,*)
x10 = cmplx(huge(real(x10))/10,huge(real(x10))/10,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(huge(real(x10))/10,1,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(huge(real(x10))/10,tiny(real(x10))*10,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(1,huge(real(x10))/10,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(tiny(real(x10))*10,huge(real(x10))/10,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(1,1,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(1,tiny(real(x10))*10,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(tiny(real(x10))*10,1,kind(x10))
write(*,*) x10, cabs(x10)
x10 = cmplx(tiny(real(x10))*10,tiny(real(x10))*10,kind(x10))
write(*,*) x10, cabs(x10)
write(*,*)
x16 = cmplx(huge(real(x16))/10,huge(real(x16))/10,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(huge(real(x16))/10,1,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(huge(real(x16))/10,tiny(real(x16))*10,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(1,huge(real(x16))/10,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(tiny(real(x16))*10,huge(real(x16))/10,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(1,1,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(1,tiny(real(x16))*10,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(tiny(real(x16))*10,1,kind(x16))
write(*,*) x16, cabs(x16)
x16 = cmplx(tiny(real(x16))*10,tiny(real(x16))*10,kind(x16))
write(*,*) x16, cabs(x16)
end program testme
C:\gfortran\clf\csqrt>gfortran cabs_test.f90 -ocabs_test
C:\gfortran\clf\csqrt>cabs_test
( 3.40282347E+37, 3.40282347E+37) 4.81231910E+37
( 3.40282347E+37, 1.0000000 ) 3.40282347E+37
( 3.40282347E+37, 1.17549435E-37) 3.40282347E+37
( 1.0000000 , 3.40282347E+37) 3.40282347E+37
( 1.17549435E-37, 3.40282347E+37) 3.40282347E+37
( 1.0000000 , 1.0000000 ) 1.4142135
( 1.0000000 , 1.17549435E-37) 1.0000000
( 1.17549435E-37, 1.0000000 ) 1.0000000
( 1.17549435E-37, 1.17549435E-37) 1.66240005E-37
( 1.79769313486231578E+307, 1.79769313486231578E+307)
2.54232201230729257E+307
( 1.79769313486231578E+307, 1.0000000000000000 )
1.79769313486231578E+307
( 1.79769313486231578E+307, 2.22507385850720138E-307)
1.79769313486231578E+307
( 1.0000000000000000 , 1.79769313486231578E+307)
1.79769313486231578E+307
( 2.22507385850720138E-307, 1.79769313486231578E+307)
1.79769313486231578E+307
( 1.0000000000000000 , 1.0000000000000000 ) 1.4142135623730951
( 1.0000000000000000 , 2.22507385850720138E-307) 1.0000000000000000
( 2.22507385850720138E-307, 1.0000000000000000 ) 1.0000000000000000
( 2.22507385850720138E-307, 2.22507385850720138E-307)
3.14672962798271743E-307
( n.inite 85850720138262E+0000, n.inite 85850720138263E+0000) n.inite OD
?
E+0000
( n.inite 85850720138264E+0000, 1.0000000000000000000 ) n.inite OD
?
E+0000
( n.inite 00000000000000E+0000, 0.0000000000000-(00000E-4932) n.inite OD
?
E+0000
( 1.0000000000000000000 , n.inite 00000000000000E+0000) n.inite OD
?
E+0000
( 0.0000000000000-(00000E-4932, n.inite 0000000-(00000E-4932) n.inite OD
?
E+0000
( 1.0000000000000000000 , 1.0000000000000000000 )
1.4142135623730
950487
( 1.0000000000000000000 , 0.0000000000000-(00000E-4932)
1.0000000000000
000000
( 0.0000000000000-(00000E-4932, 1.0000000000000000000 )
1.0000000000000
000000
( 0.0000000000000-(00000E-4932, 0.0000000000000-(00000E-4932)
0.0000000000000-
(00000E-4932
( 1.18973149535723176508575932662800704E+4931,
1.189731495357231765085759326628
00704E+4931) 1.68253441631662012728269141999214698E+4931
( 1.18973149535723176508575932662800704E+4931,
1.00000000000000000000000000000
00000 ) 1.18973149535723176508575932662800704E+4931
( 1.18973149535723176508575932662800704E+4931,
3.362103143112093506262677817321
75260E-4931) 1.18973149535723176508575932662800704E+4931
( 1.0000000000000000000000000000000000 ,
1.189731495357231765085759326628
00704E+4931) 1.18973149535723176508575932662800704E+4931
( 3.36210314311209350626267781732175260E-4931,
1.189731495357231765085759326628
00704E+4931) 1.18973149535723176508575932662800704E+4931
( 1.0000000000000000000000000000000000 ,
1.00000000000000000000000000000
00000 ) 1.4142135623730950488016887242096982
( 1.0000000000000000000000000000000000 ,
3.362103143112093506262677817321
75260E-4931) 1.0000000000000000000000000000000000
( 3.36210314311209350626267781732175260E-4931,
1.00000000000000000000000000000
00000 ) 1.0000000000000000000000000000000000
( 3.36210314311209350626267781732175260E-4931,
3.362103143112093506262677817321
75260E-4931) 4.75473186308633355902452847925549340E-4931
So no overflow or underflow, just some bugs in real(kind=10)
intrinsics. Now let's apply this to mandelbrot:
Place module m8 near the top of mandelbrot.f90:
! gfortran -g mandelbrot.f90 `pkg-config --cflags --libs gtk+-2.0`
! Contributed by Jerry DeLisle and Vincent Magnin
module m8
use ISO_FORTRAN_ENV
implicit none
private
integer, parameter :: wp = REAL_KINDS(2)
public cabs
interface cabs
module procedure cabs_gen
end interface cabs
contains
include 'cabs.f90'
end module m8
USE it in program mandelbrot:
program mandelbrot
use iso_c_binding, only: c_ptr, c_null_char, c_null_ptr, c_funloc
use m8
use handlers
And switch to our external function:
! do while ((k <= itermax) .and. (abs(z)<2d0))
do while ((k <= itermax) .and. (cabs(z)<2d0))
! do while ((k <= itermax) .and. (real(conjg(z)*z)<4d0))
Compile, link, and run. My result was that execution time
went down from 3 minutes to 30 seconds. Forgot important
modification of cabs.f90:
C:\gfortran\GTK\test1>type cabs.f90
! File: cabs.f90
! Public domain 2011 James Van Buskirk
elemental function cabs_gen(z)
implicit none
complex(wp),intent(in),value :: z
real(wp),parameter :: x = 0
real(wp),parameter :: tolhi = nearest(sqrt(huge(x)/2),-1.0)
real(wp),parameter :: tollo = nearest(sqrt(tiny(x)/epsilon(x)),1.0)
real(wp),parameter :: mulhi = scale(1.0_wp,exponent(tolhi/huge(x))-1)
real(wp),parameter :: divhi = scale(1.0_wp,1-exponent(tolhi/huge(x)))
real(wp),parameter :: mullo = scale(1.0_wp,exponent(tollo/tiny(x))+1)
real(wp),parameter :: divlo = scale(1.0_wp,-1-exponent(tollo/tiny(x)))
real(wp) cabs_gen
if(abs(real(z)) > tolhi) then
cabs_gen = divhi*sqrt((mulhi*real(z))**2+(mulhi*aimag(z))**2)
else if(abs(real(z)) < tollo) then
if(abs(aimag(z)) < tollo) then
cabs_gen = divlo*sqrt((mullo*real(z))**2+(mullo*aimag(z))**2)
else if(abs(aimag(z)) > tolhi) then
cabs_gen = abs(aimag(z))
else
cabs_gen = sqrt(real(z)**2+aimag(z)**2)
end if
else
if(abs(aimag(z)) > tolhi) then
cabs_gen = divhi*sqrt((mulhi*real(z))**2+(mulhi*aimag(z))**2)
else
cabs_gen = sqrt(real(z)**2+aimag(z)**2)
end if
end if
end function cabs_gen
Vincent
made some tests with my mandelbrot program (without GUI) using
gfortran mandelbrot.f90 (version 4.6.0 under Linux) and obtained the
following results:
do while ((k <= itermax) .and. (abs(z)<2d0))
41.3 s
do while ((k <= itermax) .and. (real(conjg(z)*z)<4d0))
10.3 s
do while ((k <= itermax) .and. ((real(z)*real(z)
+aimag(z)*aimag(z))<4d0))
10.1 s
do while ((k <= itermax) .and. ((real(z)**2+aimag(z)**2)<4d0))
9.72 s
and also:
do while ((k <= itermax) .and. (sqrt(real(z)**2+aimag(z)**2)<2d0))
14.2 s