The error I get when I try to use a global reference to the array is now: Assertion failed in file commutil.c at line 252: comm->node_comm == ((void *)0) internal ABORT - process 0
I have spent some time in the mailing list archives looking for examples of how to write a material function such that an array of dielectric values can be used. I did successfully create a hdf5 file with the corresponding values using matlab. I have added links to a number of related threads in the mailing list and appended version of holey-wvg-cavity.cpp that reads dielectric values from eps.txt (available from http://www.doe.carleton.ca/~kmedri/research/eps.txt) successfully into an array. It also includes code to return the corresponding value when provided the vector and the array. I conclude from this that hdf5 route is not the way to go: "Reloading eps files" http://thread.gmane.org/gmane.comp.science.electromagnetism.meep.general/964 /focus=1159 I would appreciate a code sample from: "What is the most efficient way to import data to meep?" http://thread.gmane.org/gmane.comp.science.electromagnetism.meep.general/228 4/focus=2286 as well as from: "building up non periodic structure using c++ interface" http://thread.gmane.org/gmane.comp.science.electromagnetism.meep.general/296 1 I am interested in writing a material function: "Arbitrary Shapes in Meep" http://thread.gmane.org/gmane.comp.science.electromagnetism.meep.general/131 /focus=132 This is something to work towards: "Voxel file interface for Meep" http://thread.gmane.org/gmane.comp.science.electromagnetism.meep.general/213 8/focus=2146 #include "meep.hpp" using namespace meep; #include <fstream> #include <iostream> /* Based on holey-wvg-cavity.ctl, etc. from the Meep Tutorial */ //Some parameters to describe the geometry: const double epswg = 13; //dielectric constant of waveguide const double w = 1.2; //width of waveguide const double r = 0.36; //radius of holes const double d = 1.4; //defect spacing const double l=1; // (ordinary spacing = 1) const double N = 3; //number of holes on either side of defect //N is not used const double sy = 6; //size of cell in y direction (perpendicular to wvg.) const double pad = 2; //padding between last hole and PML edge const double dpml = 1; //PML thickness const double sx = 2*(pad+dpml+N) + d - 1; //size of cell in x direction double darray[239][495]; //double diearray[239][495]; complex<double> one(const vec &p) { return 1.0; } double no_holes(const vec &p) { const volume v = vol2d(sx,sy,40); vec xx = p - v.center(); if (abs(xx.y()) <= w/2) return epswg; else return 1; } void pass_array(double diearray[239][495]) { master_printf("diearray[95][0]: %f\n",diearray[95][0]); master_printf("diearray[96][0]: %f\n",diearray[96][0]); } //as per 400px-Tut-holey-cavity.png double eps(const vec &p) { //once array available use this line: //return pow(diearray[round(p.y()*40)][round(p.x()*40)],2); //instead of the below code. const volume v = vol2d(sx,sy,40); vec xx = p - v.center(); vec xxN1 = vec(abs(xx.x()) - d/2,xx.y()); vec xxN2 = vec(abs(xx.x())-d/2-l,xx.y()); vec xxN3 = vec(abs(xx.x())-d/2-l*2,xx.y()); if (sqrt(pow(xxN1.x(),2)+pow(xx.y(),2))<=r) return 1; else if (sqrt(pow(xxN2.x(),2)+pow(xx.y(),2))<=r) return 1; else if (sqrt(pow(xxN3.x(),2)+pow(xx.y(),2))<=r) return 1; else if (abs(xx.y()) <= w/2) return epswg; else return 1; } void excitewg(double a, component c) { volume v = vol2d(sx,sy,a); //declaring the equivelent of odd mirror symmetry through the x axis: const symmetry S = -mirror(Y,v); structure s0(v,no_holes,pml(dpml),S); structure s(v,eps,pml(dpml),S); fields f0(&s0); fields f(&s); f0.use_real_fields(); f.use_real_fields(); //const char *dirname = make_output_directory(__FILE__); //s.set_output_directory(dirname); f0.output_hdf5(Dielectric, v.surroundings()); f.output_hdf5(Dielectric, v.surroundings()); //double fcen= 0.25; double freq= 0.25; //double df = 0.2; double fwidth = 0.2; //double fwidth = 0.2*32; double freq_min = freq - fwidth/2; double freq_max = freq + fwidth/2; int nfreq = 500; gaussian_src_time src(freq, fwidth); geometric_volume src_plane(v.center()-vec(dpml-0.5*sx,w/2), v.center()-vec(dpml-0.5*sx,-w/2)); geometric_volume flux_plane(v.center()-vec(-dpml+0.5*sx,w/2), v.center()-vec(-dpml+0.5*sx,-w/2)); f0.add_volume_source(c,src,src_plane,one,1.0); f.add_volume_source(c,src,src_plane,one,1.0); master_printf("volume sources added...\n"); dft_flux ft0 = f0.add_dft_flux_plane(flux_plane,freq_min,freq_max,nfreq); dft_flux ft = f.add_dft_flux_plane(flux_plane,freq_min,freq_max,nfreq); master_printf("simulating structure without holes...\n"); //while (f0.time() < nfreq / (freq_max-freq_min) / 2) f0.step(); while (f0.time() < 12) f0.step(); master_printf("simulating structure...\n"); //while (f.time() < nfreq / (freq_max-freq_min) / 2) f.step(); while (f.time() < 12) f.step(); double *flux = ft.flux(); double *flux0 = ft0.flux(); double *T; T = new double[nfreq]; for (int i=0; i<nfreq; ++i) T[i] = flux[i] / flux0[i]; double dfreq = (freq_max-freq_min) / (nfreq-1); FILE *TRANSMISSION_OUTPUT; TRANSMISSION_OUTPUT=master_fopen("holey-wvg-cavity.txt","wt"); master_fprintf(TRANSMISSION_OUTPUT,"tranmission:, freq, T[i]\n"); for (int i=0; i<nfreq; ++i) master_fprintf(TRANSMISSION_OUTPUT,"transmission:, %f, %f\n", freq_min+i*dfreq,T[i]); master_fclose(TRANSMISSION_OUTPUT); delete [] flux; delete [] flux0; //delete [] dirname; } int main(int argc, char **argv) { //START FILE READ char str[4000]; fstream file_op("eps.txt",ios::in); file_op >> str ; //get rid of 239 and file_op >> str ; //495 at beginning of file double dielectricarray[239][495]; for (int i=0;i<=239;i++) { for (int j=0;j<=495;j++) { file_op >> str ; dielectricarray[i][j]=strtod(str,NULL); //uncomment: //darray[i][j]=dielectricarray[i][j]; //to get //Assertion failed in file commutil.c at line 252: comm->node_comm // == ((void *)0) internal ABORT - process 0 ///mpich2/src/mpi/comm/commutil.c? //if FILE READ is after initialize mpi() and line uncommented then: //meep: v_grid_points = 60512, sum(chunks) = 0 //application called MPI_Abort(comm=0x44000000, 1) - process 0 } } file_op.close(); //END FILE READ initialize mpi(argc, argv); master_printf("holey-wvg-cavity.cpp starting.\n"); //pass dielectricarray to eps_struct(dielectricarray,vec?) pass_array(dielectricarray); excitewg(40, Ey); master_printf("holey-wvg-cavity.cpp finished.\n"); return 0; } -----Original Message----- From: meep-discuss-boun...@ab-initio.mit.edu [mailto:meep-discuss-boun...@ab-initio.mit.edu] On Behalf Of Kristian Medri Sent: February 25, 2010 12:10 PM To: meep-discuss@ab-initio.mit.edu Subject: Re: [Meep-discuss] Passing a variable to the eps funtion P.S. Creating a global variable does not cause any problems: double darray[401][401]; but if I use a for loop in main() to set the values of darray to those in dielectricarray I get the following runtime error: meep: v_grid_points = 409600, sum(chunks) = 0 application called MPI_Abort(comm=0x44000000, 1) - process 0 -----Original Message----- From: meep-discuss-boun...@ab-initio.mit.edu [mailto:meep-discuss-boun...@ab-initio.mit.edu] On Behalf Of Kristian Medri Sent: February 25, 2010 11:40 AM To: meep-discuss@ab-initio.mit.edu Subject: [Meep-discuss] Passing a variable to the eps funtion Dear Meep community, In my main{} I can call: structure s(v,eps_struct(dielectricarray),pml(dpml)); where: double eps_struct(const vec &p) I would like to pass a variable: structure s(v,eps_struct(dielectricarray),pml(dpml)); where: double eps_struct(const vec &p, double diearray[401][401]) If I do not pass anything then is finds p but when I pass dielectricarray it does not find it, any hints? Thanks! _______________________________________________ meep-discuss mailing list meep-discuss@ab-initio.mit.edu http://ab-initio.mit.edu/cgi-bin/mailman/listinfo/meep-discuss _______________________________________________ meep-discuss mailing list meep-discuss@ab-initio.mit.edu http://ab-initio.mit.edu/cgi-bin/mailman/listinfo/meep-discuss _______________________________________________ meep-discuss mailing list meep-discuss@ab-initio.mit.edu http://ab-initio.mit.edu/cgi-bin/mailman/listinfo/meep-discuss
