1 /* atlc - arbitrary transmission line calculator, for the analysis of
2 extern int number_of_workers;
3 transmission lines are directional couplers.
5 Copyright (C) 2002. Dr. David Kirkby, PhD (G8WRB).
7 This program is free software; you can redistribute it and/or
8 modify it under the terms of the GNU General Public License
9 as published by the Free Software Foundation; either package_version 2
10 of the License, or (at your option) any later package_version.
12 This program is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
17 You should have received a copy of the GNU General Public License
18 along with this program; if not, write to the Free Software
19 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
22 Dr. David Kirkby, e-mail drkirkby at gmail.com
32 #include "definitions.h"
34 extern int append_flag;
35 extern int dielectrics_to_consider_just_now, coupler;
39 extern int number_of_workers;
41 void do_fd_calculation(struct transmission_line_properties *data, size_t size, FILE *where_to_print_fp, char *inputfile_filename)
43 double capacitance_old, capacitance;
44 double velocity_of_light_in_vacuum;
47 if (data->dielectrics_in_bitmap > 1) {
48 fprintf(stderr,"\nSorry, but on the 15th October 2003 I was advised there is an\n");
49 fprintf(stderr,"error in atlc when computing systems with multiple dielectrics.\n");
50 fprintf(stderr,"So until this problem is fixed, the facility has been disabled.\n\n");
51 fprintf(stderr,"I hope to release a new version shortly without this problem.\n");
55 /* The line of best fit of non_metalic_pixels vs iterations required
56 is y=0.0011 * non_metallic_pixels + 283.
58 I'll ensure finite_difference is called about 10x by using
59 0.00011 * the number of non-metallic elements +28 as the number
60 of times finite_difference is called each time */
64 /* The following 10 lines are for a single dielectric 2 conductor line */
65 if (data->couplerQ==FALSE)
67 if(data->verbose_level >= 2)
68 printf("Solving assuming a vacuum dielectric\n");
69 capacitance=VERY_LARGE; /* Can be anything large */
70 dielectrics_to_consider_just_now=1;
71 data->dielectrics_to_consider_just_now=1;
73 do /* Start a finite calculation */
75 capacitance_old=capacitance;
77 #ifdef ENABLE_POSIX_THREADS
78 if (number_of_workers == 0)
79 capacitance=finite_difference_single_threaded();
81 capacitance=finite_difference_multi_threaded();
83 capacitance=finite_difference_single_threaded();
86 data->C_vacuum=capacitance;
88 data->L_vacuum=MU_0*EPSILON_0/capacitance; /* Same as L in *ALL* cases */
89 data->Zo_vacuum=sqrt(data->L_vacuum/data->C_vacuum); /* Standard formaul for Zo */
91 if (data->dielectrics_in_bitmap == 1) /* Just get C by simple scaling of Er */
93 data->C=capacitance*data->found_this_dielectric; /* Scaled by the single dielectric constant */
94 data->Er=data->found_this_dielectric;
98 data->Zo=sqrt(data->L_vacuum/data->C); /* Standard formula for Zo */
99 data->Zodd=sqrt(data->L_vacuum/data->C); /* Standard formula for Zo */
100 velocity_of_light_in_vacuum=1.0/(sqrt(MU_0 * EPSILON_0)); /* around 3x10^8 m/s */
101 data->velocity=1.0/pow(data->L_vacuum*data->C,0.5);
102 data->velocity_factor=data->velocity/velocity_of_light_in_vacuum;
103 data->relative_permittivity=sqrt(data->velocity_factor); /* ??? XXXXXX */
104 if(data->verbose_level > 0 ) /* Only needed if intermediate results wanted. */
105 print_data_for_two_conductor_lines(*data, where_to_print_fp, inputfile_filename);
107 } while (fabs((capacitance_old-capacitance)/capacitance_old) > data->cutoff); /* end of FD loop */
108 if(data->verbose_level >=4)
109 printf("Total of %d iterations ( %d calls to finite_difference() )\n",ITERATIONS*count,count);
111 if((data->write_binary_field_imagesQ == TRUE || data->write_bitmap_field_imagesQ == TRUE) && data->dielectrics_in_bitmap==1 )
112 write_fields_for_two_conductor_lines(inputfile_filename, *data, size);
113 if(data->verbose_level == 0 && data->dielectrics_in_bitmap==1 )
114 print_data_for_two_conductor_lines(*data, where_to_print_fp, inputfile_filename);
116 if ( data->dielectrics_in_bitmap >1)
118 /* We know the capacitance and inductance for the air spaced line
119 as we calculated it above. Howerver, whilst the inductance
120 is independant of the dielectric, the capacitance is not, so this
121 has to be recalculated, taking care not to alter the inductance
123 if(data->verbose_level >= 2)
124 printf("Now taking into account the permittivities of the different dielectrics for 2 conductors.\n");
126 dielectrics_to_consider_just_now=3; /* Any number > 1 */
127 data->dielectrics_to_consider_just_now=2; /* Any number > 1 */
129 capacitance=VERY_LARGE; /* Can be anything large */
131 do /* Start a finite calculation */
133 capacitance_old=capacitance;
134 #ifdef ENABLE_POSIX_THREADS
135 if (number_of_workers == 0)
136 capacitance=finite_difference_single_threaded();
138 capacitance=finite_difference_multi_threaded();
140 capacitance=finite_difference_single_threaded();
143 data->C_non_vacuum=capacitance;
144 data->Zo=sqrt(data->L_vacuum/data->C_non_vacuum); /* Standard formula for Zo */
145 data->velocity=1.0/pow(data->L_vacuum*data->C_non_vacuum,0.5);
146 data->velocity_factor=data->velocity/velocity_of_light_in_vacuum;
147 data->relative_permittivity=sqrt(data->velocity_factor); /* ??? XXXXXX */
148 data->Er=data->C/data->C_vacuum;
149 if(data->verbose_level > 0 ) /* Only needed if intermediate results wanted. */
150 print_data_for_two_conductor_lines(*data, where_to_print_fp, inputfile_filename);
151 } while (fabs((capacitance_old-capacitance)/capacitance_old) > data->cutoff); /* end of FD loop */
153 /* We must print the results now, but only bother if the verbose level was
154 not not incremented on the command line, otherwide there will be two duplicate
157 if (data->verbose_level == 0)
158 print_data_for_two_conductor_lines(*data, where_to_print_fp, inputfile_filename);
159 if(data->write_binary_field_imagesQ == TRUE || data->write_bitmap_field_imagesQ == TRUE)
160 write_fields_for_two_conductor_lines(inputfile_filename, *data, size);
163 else if (data->couplerQ==TRUE)
165 /* The properties of a couplers will be computed in 2 or 4 stages
166 1) Compute the odd-mode impedance, assuming a vacuum dielectric, or
167 if there is just one dielectric, that one.
169 2) Compute the odd-mode impedance, taking into account the effect of
170 multiple dielectrics, IF NECESSARY
172 at this point, the negative voltages will be turned into positive ones.
174 3) Compute the even-mode impedance, assuming a vacuum dielectric, or
175 if there is just one dielectric, that one.
177 4) Compute the even-mode impedance, taking into account the effect of
178 multiple dielectrics, IF NECESSARY */
180 /* Stage 1 - compute the odd mode impedance assuming single dielectric */
181 data->display = Z_ODD_SINGLE_DIELECTRIC;
182 dielectrics_to_consider_just_now=1;
183 data->dielectrics_to_consider_just_now=1;
185 capacitance=VERY_LARGE; /* Can be anything large */
186 if(data->verbose_level >= 2)
187 printf("Solving assuming a vacuum dielectric to compute the odd-mode impedance\n");
189 do /* Start a finite difference calculation */
191 capacitance_old=capacitance;
192 #ifdef ENABLE_POSIX_THREADS
193 if (number_of_workers == 0)
194 capacitance=finite_difference_single_threaded();
196 capacitance=finite_difference_multi_threaded();
198 capacitance=finite_difference_single_threaded();
200 data->Codd_vacuum=capacitance;
201 data->Codd=capacitance;
202 data->Lodd_vacuum=MU_0*EPSILON_0/capacitance; /* Same as L in *ALL* cases */
204 data->Zodd_vacuum=sqrt(data->Lodd_vacuum/data->Codd_vacuum); /* Standard formaul for Zodd */
206 if (data->dielectrics_in_bitmap == 1) /* Just get C by simple scaling of Er */
207 data->Codd*=data->found_this_dielectric; /* Scaled by the single dielectric constant */
210 data->Zodd=sqrt(data->Lodd_vacuum/data->Codd); /* Standard formula for Zo */
211 velocity_of_light_in_vacuum=1.0/(sqrt(MU_0 * EPSILON_0)); /* around 3x10^8 m/s */
212 /* FPE trapdata->velocity_odd=1.0/pow(data->L_vacuum*data->Codd,0.5); */
213 data->velocity_odd=1.0/pow(data->Lodd_vacuum*data->Codd,0.5);
214 data->velocity_factor_odd=data->velocity_odd/velocity_of_light_in_vacuum;
215 data->relative_permittivity_odd=sqrt(data->velocity_factor_odd); /* ??? XXXXXX */
216 data->Er_odd=data->Codd/data->Codd_vacuum;
217 data->Zdiff=2.0*data->Zodd;
218 /* Print text if uses wants it */
219 if(data->verbose_level>=1)
220 print_data_for_directional_couplers(*data, where_to_print_fp, inputfile_filename);
221 } while (fabs((capacitance_old-capacitance)/capacitance_old) > data->cutoff); /* end of FD loop */
224 mpi_receive_updated_vij_strips();
225 #endif /* ENABLE_MPI */
227 /* display bitpamps/binary files if this is the last odd-mode computation */
228 if((data->write_binary_field_imagesQ == TRUE || data->write_bitmap_field_imagesQ == TRUE) && data->dielectrics_in_bitmap==1 )
229 write_fields_for_directional_couplers(inputfile_filename, *data, size, ODD);
231 /* Stage 2 - compute the odd-mode impedance taking into account other dielectrics IF NECESSARY */
233 if ( data->dielectrics_in_bitmap >1)
235 if(data->verbose_level >= 2)
236 printf("Now taking into account the permittivities of the different dielectrics to compute Zodd.\n");
237 data->display = Z_ODD_SINGLE_DIELECTRIC;
238 capacitance=VERY_LARGE; /* Can be anything large */
240 dielectrics_to_consider_just_now=2;
241 data->dielectrics_to_consider_just_now=2;
244 mpi_send_current_data();
245 #endif /* ENABLE_MPI */
247 do /* Start a finite calculation */
249 capacitance_old=capacitance;
250 #ifdef ENABLE_POSIX_THREADS
251 if (number_of_workers == 0)
252 capacitance=finite_difference_single_threaded();
254 capacitance=finite_difference_multi_threaded();
256 capacitance=finite_difference_single_threaded();
258 data->Codd=capacitance;
259 data->Zodd=sqrt(data->Lodd_vacuum/data->Codd); /* Standard formula for Zo */
260 velocity_of_light_in_vacuum=1.0/(sqrt(MU_0 * EPSILON_0)); /* around 3x10^8 m/s */
261 data->velocity_odd=1.0/pow(data->L_vacuum*data->C,0.5);
262 data->velocity_factor_odd=data->velocity/velocity_of_light_in_vacuum;
263 data->relative_permittivity_odd=sqrt(data->velocity_factor); /* ??? XXXXXX */
264 data->Er_odd=data->Codd/data->Codd_vacuum;
265 data->Zdiff=2.0*data->Zodd;
266 if(data->verbose_level>=1)
267 print_data_for_directional_couplers(*data, where_to_print_fp, inputfile_filename);
268 } while (fabs((capacitance_old-capacitance)/capacitance_old) > data->cutoff); /* end of FD loop */
271 mpi_receive_updated_vij_strips();
272 #endif /* ENABLE_MPI */
274 if((data->write_binary_field_imagesQ == TRUE || data->write_bitmap_field_imagesQ == TRUE) && data->dielectrics_in_bitmap!=1 )
275 write_fields_for_directional_couplers(inputfile_filename, *data, size, ODD);
276 } /* end of stage 2 for couplers */
278 /* Stage 3 - compute the even-mode impedance assuming single dielectric */
280 /* Since we want the even mode impedance now, we swap all the -1V
281 metallic conductors for +1V */
283 swap_conductor_voltages();
285 data->display = Z_EVEN_SINGLE_DIELECTRIC;
286 dielectrics_to_consider_just_now=1;
287 data->dielectrics_to_consider_just_now=1;
288 if(data->verbose_level >= 2)
289 printf("Now assuming a vacuum dielectric to compute Zeven\n");
291 capacitance=VERY_LARGE; /* Can be anything large */
294 mpi_send_current_data();
295 #endif /* ENABLE_MPI */
297 do /* Start a finite difference calculation */
299 capacitance_old=capacitance;
300 #ifdef ENABLE_POSIX_THREADS
301 if (number_of_workers == 0)
302 capacitance=finite_difference_single_threaded();
304 capacitance=finite_difference_multi_threaded();
306 capacitance=finite_difference_single_threaded();
309 data->Ceven_vacuum=capacitance;
310 data->Ceven=capacitance;
311 data->Leven_vacuum=MU_0*EPSILON_0/capacitance; /* Same as L in *ALL* cases */
313 data->Zeven_vacuum=sqrt(data->Leven_vacuum/data->Ceven_vacuum); /* Standard formaul for Zodd */
315 if (data->dielectrics_in_bitmap == 1) /* Just get C by simple scaling of Er */
316 data->Ceven*=data->found_this_dielectric; /* Scaled by the single dielectric constant */
319 data->Zeven=sqrt(data->Leven_vacuum/data->Ceven); /* Standard formula for Zo */
320 velocity_of_light_in_vacuum=1.0/(sqrt(MU_0 * EPSILON_0)); /* around 3x10^8 m/s */
321 data->velocity_even=1.0/pow(data->Leven_vacuum*data->Ceven,0.5);
322 data->velocity_factor_even=data->velocity_even/velocity_of_light_in_vacuum;
323 data->relative_permittivity_even=sqrt(data->velocity_factor_even); /* ??? XXXXXX */
324 data->Er_even=data->Ceven/data->Ceven_vacuum;
325 data->Zcomm=data->Zeven/2.0;
326 data->Zo=sqrt(data->Zodd * data->Zeven);
327 if(data->verbose_level>=1)
328 print_data_for_directional_couplers(*data, where_to_print_fp, inputfile_filename);
329 /* display bitpamps/binary files if this is the last even-mode computation */
330 } while (fabs((capacitance_old-capacitance)/capacitance_old) > data->cutoff); /* end of FD loop */
332 if((data->write_binary_field_imagesQ == TRUE || data->write_bitmap_field_imagesQ == TRUE) && data->dielectrics_in_bitmap==1)
333 write_fields_for_directional_couplers(inputfile_filename, *data, size, EVEN);
335 capacitance=VERY_LARGE; /* Can be anything large */
336 /* Stage 4 - compute the even-mode impedance assuming multiple dielectrics IF NECESSARY */
337 if ( data->dielectrics_in_bitmap >1)
339 dielectrics_to_consider_just_now=2;
340 data->dielectrics_to_consider_just_now=2;
341 if(data->verbose_level >= 2)
342 printf("Now taking into account the permittivities of the different dielectrics to compute Zeven\n");
345 mpi_send_current_data();
346 #endif /* ENABLE_MPI */
348 do /* Start a finite calculation */
350 capacitance_old=capacitance;
351 #ifdef ENABLE_POSIX_THREADS
352 if (number_of_workers == 0)
353 capacitance=finite_difference_single_threaded();
355 capacitance=finite_difference_multi_threaded();
357 capacitance=finite_difference_single_threaded();
359 data->Ceven=capacitance;
360 data->Zeven=sqrt(data->Leven_vacuum/data->Ceven); /* Standard formula for Zo */
361 velocity_of_light_in_vacuum=1.0/(sqrt(MU_0 * EPSILON_0)); /* around 3x10^8 m/s */
362 data->velocity_even=1.0/pow(data->L_vacuum*data->C,0.5);
363 data->velocity_factor_even=data->velocity/velocity_of_light_in_vacuum;
364 data->relative_permittivity_even=sqrt(data->velocity_factor); /* ??? XXXXXX */
365 data->Er_even=data->Ceven/data->Ceven_vacuum;
366 data->Zdiff=2.0*data->Zodd;
367 data->Zcomm=data->Zeven/2.0;
368 data->Zo=sqrt(data->Zeven*data->Zodd);
369 if(data->verbose_level>=1)
370 print_data_for_directional_couplers(*data, where_to_print_fp, inputfile_filename);
371 } while (fabs((capacitance_old-capacitance)/capacitance_old) > data->cutoff); /* end of FD loop */
373 if(data->write_binary_field_imagesQ == TRUE || data->write_bitmap_field_imagesQ == TRUE)
374 write_fields_for_directional_couplers(inputfile_filename, *data, size, EVEN);
375 } /* end of stage 4 */
376 /* Print the results if the verbose level was 0 (no -v flag(s) ). */
377 if (data->verbose_level == 0)
379 /* We need to print the data. The next function will only print if
380 the verbose_level is 1 or more, so I'll fix it at one. Then we print
381 the final results and exit. */
382 data->verbose_level=1;
383 data->display = Z_EVEN_SINGLE_DIELECTRIC;
384 print_data_for_directional_couplers(*data, where_to_print_fp, inputfile_filename);
386 } /* end of if couplers */