/[escript]/trunk/finley/src/Mesh_createNodeFileMappings.c
ViewVC logotype

Contents of /trunk/finley/src/Mesh_createNodeFileMappings.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1735 - (show annotations)
Fri Aug 29 00:18:26 2008 UTC (10 years, 7 months ago) by gross
File MIME type: text/plain
File size: 18230 byte(s)
This fixes a bug in the callculation of communication tables in
particular for the data to be send away. The problem was that the
previously used method based on the idea to use the elements stored on a
processor to decide what to send away is not suitable for the reduced
DOF. Now each processor sends a list of the DOFs it is expecting to its
naighbours. 



1
2 /* $Id:$ */
3
4 /*******************************************************
5 *
6 * Copyright 2003-2007 by ACceSS MNRF
7 * Copyright 2007 by University of Queensland
8 *
9 * http://esscc.uq.edu.au
10 * Primary Business: Queensland, Australia
11 * Licensed under the Open Software License version 3.0
12 * http://www.opensource.org/licenses/osl-3.0.php
13 *
14 *******************************************************/
15
16 /**************************************************************/
17
18 /* Finley: NodeFile : creates the mappings using the indexReducedNodes */
19 /* no distribution is happening */
20
21 /**************************************************************/
22
23 #include "Mesh.h"
24 #define UNUSED -1
25
26 #define BOUNDS_CHECK 1
27
28 /**************************************************************/
29
30 void Mesh_createDOFMappingAndCoupling(Finley_Mesh* in, bool_t use_reduced_elements)
31 {
32 index_t min_DOF, max_DOF, *shared=NULL, *offsetInShared=NULL, *locDOFMask=NULL, i, k, myFirstDOF, myLastDOF, *nodeMask=NULL, firstDOF, lastDOF, *globalDOFIndex, *wanted_DOFs=NULL;
33 dim_t mpiSize, len_loc_dof, numNeighbors, n, lastn, numNodes,*rcv_len=NULL, *snd_len=NULL, count;
34 Paso_MPI_rank myRank,p,p_min,p_max, *neighbor=NULL;
35 Paso_SharedComponents *rcv_shcomp=NULL, *snd_shcomp=NULL;
36 Finley_NodeMapping *this_mapping=NULL;
37 Paso_Connector* this_connector=NULL;
38 Paso_Distribution* dof_distribution;
39 Paso_MPIInfo *mpi_info = in->MPIInfo;
40 #ifdef PASO_MPI
41 MPI_Request* mpi_requests=NULL;
42 MPI_Status* mpi_stati=NULL;
43 #else
44 int *mpi_requests=NULL, *mpi_stati=NULL;
45 #endif
46
47 numNodes=in->Nodes->numNodes;
48 if (use_reduced_elements) {
49 dof_distribution=in->Nodes->reducedDegreesOfFreedomDistribution;
50 globalDOFIndex=in->Nodes->globalReducedDOFIndex;
51 } else {
52 dof_distribution=in->Nodes->degreesOfFreedomDistribution;
53 globalDOFIndex=in->Nodes->globalDegreesOfFreedom;
54 }
55
56 mpiSize=mpi_info->size;
57 myRank=mpi_info->rank;
58
59 min_DOF=Finley_Util_getFlaggedMinInt(1,numNodes,globalDOFIndex,-1);
60 max_DOF=Finley_Util_getFlaggedMaxInt(1,numNodes,globalDOFIndex,-1);
61
62 p_min=mpiSize;
63 p_max=-1;
64
65 for (p=0; p<mpiSize; ++p) {
66 if (dof_distribution->first_component[p]<=min_DOF) p_min=p;
67 if (dof_distribution->first_component[p]<=max_DOF) p_max=p;
68 }
69
70 len_loc_dof=max_DOF-min_DOF+1;
71 myFirstDOF=Paso_Distribution_getFirstComponent(dof_distribution);
72 myLastDOF=Paso_Distribution_getLastComponent(dof_distribution);
73 if (! ((min_DOF<=myFirstDOF) && (myLastDOF-1<=max_DOF)) ) {
74 Finley_setError(SYSTEM_ERROR,"Local elements do not span local degrees of freedom.");
75 return;
76 }
77 rcv_len=TMPMEMALLOC(mpiSize,dim_t);
78 snd_len=TMPMEMALLOC(mpiSize,dim_t);
79 #ifdef PASO_MPI
80 mpi_requests=MEMALLOC(mpiSize*2,MPI_Request);
81 mpi_stati=MEMALLOC(mpiSize*2,MPI_Status);
82 #else
83 mpi_requests=MEMALLOC(mpiSize*2,int);
84 mpi_stati=MEMALLOC(mpiSize*2,int);
85 #endif
86 wanted_DOFs=TMPMEMALLOC(numNodes,index_t);
87 nodeMask=TMPMEMALLOC(numNodes,index_t);
88 neighbor=TMPMEMALLOC(mpiSize,Paso_MPI_rank);
89 shared=TMPMEMALLOC(numNodes*(p_max-p_min+1),index_t);
90 offsetInShared=TMPMEMALLOC(mpiSize+1,index_t);
91 locDOFMask=TMPMEMALLOC(len_loc_dof, index_t);
92 if (! ( Finley_checkPtr(neighbor) || Finley_checkPtr(shared) || Finley_checkPtr(offsetInShared) || Finley_checkPtr(locDOFMask) ||
93 Finley_checkPtr(nodeMask) || Finley_checkPtr(rcv_len) || Finley_checkPtr(snd_len) || Finley_checkPtr(mpi_requests) || Finley_checkPtr(mpi_stati) ||
94 Finley_checkPtr(mpi_stati) )) {
95
96 memset(rcv_len,0,sizeof(dim_t)*mpiSize);
97 #pragma omp parallel
98 {
99 #pragma omp for private(i) schedule(static)
100 for (i=0;i<len_loc_dof;++i) locDOFMask[i]=UNUSED;
101 #pragma omp for private(i) schedule(static)
102 for (i=0;i<numNodes;++i) nodeMask[i]=UNUSED;
103 #pragma omp for private(i,k) schedule(static)
104 for (i=0;i<numNodes;++i) {
105 k=globalDOFIndex[i];
106 if (k>-1) {
107 locDOFMask[k-min_DOF]=UNUSED-1;
108 #ifdef BOUNDS_CHECK
109 if ((k-min_DOF) >= len_loc_dof) { printf("BOUNDS_CHECK %s %d i=%d k=%d min_DOF=%d\n", __FILE__, __LINE__, i, k, min_DOF); exit(1); }
110 #endif
111 }
112 }
113
114 #pragma omp for private(i) schedule(static)
115 for (i=myFirstDOF-min_DOF;i<myLastDOF-min_DOF;++i) {
116 locDOFMask[i]=i-myFirstDOF+min_DOF;
117 #ifdef BOUNDS_CHECK
118 if (i < 0 || i >= len_loc_dof) { printf("BOUNDS_CHECK %s %d i=%d\n", __FILE__, __LINE__, i); exit(1); }
119 #endif
120 }
121 }
122
123 numNeighbors=0;
124 n=0;
125 lastn=n;
126 for (p=p_min;p<=p_max;++p) {
127 firstDOF=MAX(min_DOF,dof_distribution->first_component[p]);
128 lastDOF=MIN(max_DOF+1,dof_distribution->first_component[p+1]);
129 if (p != myRank) {
130 for (i=firstDOF-min_DOF;i<lastDOF-min_DOF;++i) {
131 #ifdef BOUNDS_CHECK
132 if (i < 0 || i >= len_loc_dof) { printf("BOUNDS_CHECK %s %d p=%d i=%d\n", __FILE__, __LINE__, p, i); exit(1); }
133 #endif
134 if (locDOFMask[i] == UNUSED-1) {
135 locDOFMask[i]=myLastDOF-myFirstDOF+n;
136 wanted_DOFs[n]=i+min_DOF;
137 ++n;
138 }
139 }
140 if (n>lastn) {
141 rcv_len[p]=n-lastn;
142 neighbor[numNeighbors]=p;
143 #ifdef BOUNDS_CHECK
144 if (numNeighbors < 0 || numNeighbors >= mpiSize+1) { printf("BOUNDS_CHECK %s %d p=%d numNeighbors=%d n=%d\n", __FILE__, __LINE__, p, numNeighbors, n); exit(1); }
145 #endif
146 offsetInShared[numNeighbors]=lastn;
147 numNeighbors++;
148 lastn=n;
149 }
150 }
151 }
152 #ifdef BOUNDS_CHECK
153 if (numNeighbors < 0 || numNeighbors >= mpiSize+1) { printf("BOUNDS_CHECK %s %d numNeighbors=%d\n", __FILE__, __LINE__, numNeighbors); exit(1); }
154 #endif
155 offsetInShared[numNeighbors]=lastn;
156
157 /* assign new DOF labels to nodes */
158 #pragma omp parallel for private(i) schedule(static)
159 for (i=0;i<numNodes;++i) {
160 k=globalDOFIndex[i];
161 if (k>-1) nodeMask[i]=locDOFMask[k-min_DOF];
162 }
163
164 /* now we can set the mapping from nodes to local DOFs */
165 this_mapping=Finley_NodeMapping_alloc(numNodes,nodeMask,UNUSED);
166 /* define how to get DOF values for controlled bu other processors */
167 #ifdef BOUNDS_CHECK
168 for (i=0;i<offsetInShared[numNeighbors];++i) {
169 if (i < 0 || i >= numNodes*(p_max-p_min+1)) { printf("BOUNDS_CHECK %s %d i=%d\n", __FILE__, __LINE__, i); exit(1); }
170 }
171 #endif
172 #pragma omp parallel for private(i) schedule(static)
173 for (i=0;i<offsetInShared[numNeighbors];++i) shared[i]=myLastDOF-myFirstDOF+i;
174
175 rcv_shcomp=Paso_SharedComponents_alloc(myLastDOF-myFirstDOF,numNeighbors,neighbor,shared,offsetInShared,1,0,mpi_info);
176
177 /*
178 * now we build the sender
179 */
180 #ifdef PASO_MPI
181 MPI_Alltoall(rcv_len,1,MPI_INT,snd_len,1,MPI_INT,mpi_info->comm);
182 #else
183 for (p=0;p<mpiSize;++p) snd_len[p]=rcv_count[p];
184 #endif
185 count=0;
186 for (p=0;p<rcv_shcomp->numNeighbors;p++) {
187 #ifdef PASO_MPI
188 MPI_Isend(&(wanted_DOFs[rcv_shcomp->offsetInShared[p]]), rcv_shcomp->offsetInShared[p+1]-rcv_shcomp->offsetInShared[p],
189 MPI_INT,rcv_shcomp->neighbor[p],mpi_info->msg_tag_counter+myRank,mpi_info->comm,&mpi_requests[count]);
190 #endif
191 count++;
192 }
193 n=0;
194 numNeighbors=0;
195 for (p=0;p<mpiSize;p++) {
196 if (snd_len[p] > 0) {
197 #ifdef PASO_MPI
198 MPI_Irecv(&(shared[n]),snd_len[p],
199 MPI_INT, p, mpi_info->msg_tag_counter+p, mpi_info->comm, &mpi_requests[count]);
200 #endif
201 count++;
202 neighbor[numNeighbors]=p;
203 offsetInShared[numNeighbors]=n;
204 numNeighbors++;
205 n+=snd_len[p];
206 }
207 }
208 mpi_info->msg_tag_counter+=mpi_info->size;
209 offsetInShared[numNeighbors]=n;
210 #ifdef PASO_MPI
211 MPI_Waitall(count,mpi_requests,mpi_stati);
212 #endif
213 /* map global ids to local id's */
214 #pragma omp parallel for private(i) schedule(static)
215 for (i=0;i<offsetInShared[numNeighbors];++i) {
216 shared[i]=locDOFMask[shared[i]-min_DOF];
217 }
218
219 snd_shcomp=Paso_SharedComponents_alloc(myLastDOF-myFirstDOF,numNeighbors,neighbor,shared,offsetInShared,1,0,dof_distribution->mpi_info);
220
221 if (Finley_noError()) this_connector=Paso_Connector_alloc(snd_shcomp,rcv_shcomp);
222 /* assign new DOF labels to nodes */
223 Paso_SharedComponents_free(rcv_shcomp);
224 Paso_SharedComponents_free(snd_shcomp);
225 }
226 TMPMEMFREE(rcv_len);
227 TMPMEMFREE(snd_len);
228 TMPMEMFREE(mpi_requests);
229 TMPMEMFREE(mpi_stati);
230 TMPMEMFREE(wanted_DOFs);
231 TMPMEMFREE(nodeMask);
232 TMPMEMFREE(neighbor);
233 TMPMEMFREE(shared);
234 TMPMEMFREE(offsetInShared);
235 TMPMEMFREE(locDOFMask);
236 if (Finley_noError()) {
237 if (use_reduced_elements) {
238 in->Nodes->reducedDegreesOfFreedomMapping=this_mapping;
239 in->Nodes->reducedDegreesOfFreedomConnector=this_connector;
240 } else {
241 in->Nodes->degreesOfFreedomMapping=this_mapping;
242 in->Nodes->degreesOfFreedomConnector=this_connector;
243 }
244 } else {
245 Finley_NodeMapping_free(this_mapping);
246 Paso_Connector_free(this_connector);
247
248 }
249 }
250
251 void Finley_Mesh_createMappings(Finley_Mesh* mesh, index_t* distribution) {
252 int i;
253 index_t *maskReducedNodes=NULL, *indexReducedNodes=NULL;
254 dim_t numReducedNodes;
255
256 maskReducedNodes=TMPMEMALLOC(mesh->Nodes->numNodes,index_t);
257 indexReducedNodes=TMPMEMALLOC(mesh->Nodes->numNodes,index_t);
258
259 if (! ( Finley_checkPtr(maskReducedNodes) || Finley_checkPtr(indexReducedNodes) ) ) {
260 #pragma omp parallel for private(i) schedule(static)
261 for (i=0;i<mesh->Nodes->numNodes;++i) maskReducedNodes[i]=-1;
262 Finley_Mesh_markNodes(maskReducedNodes,0,mesh,TRUE);
263 numReducedNodes=Finley_Util_packMask(mesh->Nodes->numNodes,maskReducedNodes,indexReducedNodes);
264 if (Finley_noError()) Finley_Mesh_createNodeFileMappings(mesh,numReducedNodes,indexReducedNodes,distribution);
265 }
266
267 TMPMEMFREE(maskReducedNodes);
268 TMPMEMFREE(indexReducedNodes);
269 }
270
271 void Finley_Mesh_createNodeFileMappings(Finley_Mesh* in, dim_t numReducedNodes, index_t* indexReducedNodes, index_t* dof_first_component) {
272
273
274 index_t myFirstDOF, myLastDOF, myFirstNode, myLastNode, *reduced_dof_first_component=NULL, *nodeMask=NULL,
275 *reduced_nodes_first_component=NULL, *nodes_first_component=NULL,k,
276 *maskMyReducedDOF=NULL, *indexMyReducedDOF=NULL, *maskMyReducedNodes=NULL, *indexMyReducedNodes=NULL;
277 dim_t myNumDOF, myNumNodes, myNumReducedNodes, myNumReducedDOF, globalNumReducedNodes, globalNumReducedDOF,i,mpiSize, globalNumNodes;
278 Paso_MPI_rank myRank;
279
280 mpiSize=in->Nodes->MPIInfo->size;
281 myRank=in->Nodes->MPIInfo->rank;
282 /* mark the nodes used by the reduced mesh */
283
284 reduced_dof_first_component=TMPMEMALLOC(mpiSize+1,index_t);
285 reduced_nodes_first_component=TMPMEMALLOC(mpiSize+1,index_t);
286 nodes_first_component=TMPMEMALLOC(mpiSize+1,index_t);
287
288 if (! ( Finley_checkPtr(reduced_dof_first_component) || Finley_checkPtr(reduced_nodes_first_component) || Finley_checkPtr(nodes_first_component) ) ) {
289
290 globalNumNodes=Finley_NodeFile_maxGlobalNodeIDIndex(in->Nodes)+1;
291 Paso_MPIInfo_setDistribution(in->Nodes->MPIInfo,0,globalNumNodes-1,nodes_first_component);
292
293 myFirstDOF=dof_first_component[myRank];
294 myLastDOF=dof_first_component[myRank+1];
295 myNumDOF=myLastDOF-myFirstDOF;
296 myFirstNode=nodes_first_component[myRank];
297 myLastNode=nodes_first_component[myRank+1];
298 myNumNodes=myLastNode-myFirstNode;
299
300 maskMyReducedDOF=TMPMEMALLOC(myNumDOF,index_t);
301 indexMyReducedDOF=TMPMEMALLOC(myNumDOF,index_t);
302 maskMyReducedNodes=TMPMEMALLOC(myNumNodes,index_t);
303 indexMyReducedNodes=TMPMEMALLOC(myNumNodes,index_t);
304
305 if (! ( Finley_checkPtr(maskMyReducedDOF) || Finley_checkPtr(indexMyReducedDOF) || Finley_checkPtr(maskMyReducedNodes) || Finley_checkPtr(indexMyReducedNodes) ) ) {
306
307 #pragma omp parallel private(i)
308 {
309 #pragma omp for schedule(static)
310 for (i=0;i<myNumNodes;++i) maskMyReducedNodes[i]=-1;
311 #pragma omp for schedule(static)
312 for (i=0;i<myNumDOF;++i) maskMyReducedDOF[i]=-1;
313 #pragma omp for private(k) schedule(static)
314 for (i=0;i<numReducedNodes;++i) {
315 k=in->Nodes->globalNodesIndex[indexReducedNodes[i]];
316 if ( (k>=myFirstNode) && (myLastNode>k) ) maskMyReducedNodes[k-myFirstNode]=i;
317 k=in->Nodes->globalDegreesOfFreedom[indexReducedNodes[i]];
318 if ( (k>=myFirstDOF) && (myLastDOF>k) ) {
319 maskMyReducedDOF[k-myFirstDOF]=i;
320 }
321 }
322 }
323 myNumReducedNodes=Finley_Util_packMask(myNumNodes,maskMyReducedNodes,indexMyReducedNodes);
324 myNumReducedDOF=Finley_Util_packMask(myNumDOF,maskMyReducedDOF,indexMyReducedDOF);
325
326 #ifdef PASO_MPI
327 MPI_Allgather(&myNumReducedNodes,1,MPI_INT,reduced_nodes_first_component,1,MPI_INT,in->Nodes->MPIInfo->comm);
328 MPI_Allgather(&myNumReducedDOF,1,MPI_INT,reduced_dof_first_component,1,MPI_INT,in->Nodes->MPIInfo->comm);
329 #else
330 reduced_nodes_first_component[0]=myNumReducedNodes;
331 reduced_dof_first_component[0]=myNumReducedDOF;
332 #endif
333 globalNumReducedNodes=0;
334 globalNumReducedDOF=0;
335 for (i=0;i<mpiSize;++i) {
336 k=reduced_nodes_first_component[i];
337 reduced_nodes_first_component[i]=globalNumReducedNodes;
338 globalNumReducedNodes+=k;
339
340 k=reduced_dof_first_component[i];
341 reduced_dof_first_component[i]=globalNumReducedDOF;
342 globalNumReducedDOF+=k;
343 }
344 reduced_nodes_first_component[mpiSize]=globalNumReducedNodes;
345 reduced_dof_first_component[mpiSize]=globalNumReducedDOF;
346 /* ==== distribution of Nodes ===============================*/
347 Paso_MPIInfo_setDistribution(in->Nodes->MPIInfo,0,globalNumNodes-1,nodes_first_component);
348 in->Nodes->nodesDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,nodes_first_component,1,0);
349
350 /* ==== distribution of Nodes ===============================*/
351 in->Nodes->degreesOfFreedomDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,dof_first_component,1,0);
352
353 /* ==== distribution of reduced Nodes ===============================*/
354 reduced_nodes_first_component[mpiSize]=globalNumReducedNodes;
355 in->Nodes->reducedNodesDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,reduced_nodes_first_component,1,0);
356
357 /* ==== distribution of reduced DOF ===============================*/
358 in->Nodes->reducedDegreesOfFreedomDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,reduced_dof_first_component,1,0);
359 }
360 TMPMEMFREE(maskMyReducedDOF);
361 TMPMEMFREE(indexMyReducedDOF);
362 TMPMEMFREE(maskMyReducedNodes);
363 TMPMEMFREE(indexMyReducedNodes);
364 }
365 TMPMEMFREE(reduced_dof_first_component);
366 TMPMEMFREE(reduced_nodes_first_component);
367 TMPMEMFREE(nodes_first_component);
368
369 nodeMask=TMPMEMALLOC(in->Nodes->numNodes,index_t);
370 if (! Finley_checkPtr(nodeMask) && Finley_noError()) {
371
372 /* ==== nodes mapping which is a dummy structure ======== */
373 #pragma omp parallel for private(i) schedule(static)
374 for (i=0;i<in->Nodes->numNodes;++i) nodeMask[i]=i;
375 in->Nodes->nodesMapping=Finley_NodeMapping_alloc(in->Nodes->numNodes,nodeMask,UNUSED);
376
377 /* ==== mapping between nodes and reduced nodes ========== */
378 #pragma omp parallel for private(i) schedule(static)
379 for (i=0;i<in->Nodes->numNodes;++i) nodeMask[i]=UNUSED;
380 #pragma omp parallel for private(i) schedule(static)
381 for (i=0;i<numReducedNodes;++i) nodeMask[indexReducedNodes[i]]=i;
382 in->Nodes->reducedNodesMapping=Finley_NodeMapping_alloc(in->Nodes->numNodes,nodeMask,UNUSED);
383
384 }
385 TMPMEMFREE(nodeMask);
386 /* ==== mapping between nodes and DOFs + DOF connector ========== */
387 if ( Finley_noError()) Mesh_createDOFMappingAndCoupling(in,FALSE);
388 /* ==== mapping between nodes and reduced DOFs + reduced DOF connector ========== */
389 if ( Finley_noError()) Mesh_createDOFMappingAndCoupling(in,TRUE);
390
391 /* get the Ids for DOFs and reduced nodes */
392 if (Finley_noError()) {
393 #pragma omp parallel private(i)
394 {
395 #pragma omp for
396 for (i=0;i<in->Nodes->reducedNodesMapping->numTargets;++i) in->Nodes->reducedNodesId[i]=in->Nodes->Id[in->Nodes->reducedNodesMapping->map[i]];
397 #pragma omp for
398 for (i=0;i<in->Nodes->degreesOfFreedomMapping->numTargets;++i) in->Nodes->degreesOfFreedomId[i]=in->Nodes->Id[in->Nodes->degreesOfFreedomMapping->map[i]];
399 #pragma omp for
400 for (i=0;i<in->Nodes->reducedDegreesOfFreedomMapping->numTargets;++i) in->Nodes->reducedDegreesOfFreedomId[i]=in->Nodes->Id[in->Nodes->reducedDegreesOfFreedomMapping->map[i]];
401 }
402 } else {
403 Finley_NodeMapping_free(in->Nodes->nodesMapping);
404 Finley_NodeMapping_free(in->Nodes->reducedNodesMapping);
405 Finley_NodeMapping_free(in->Nodes->degreesOfFreedomMapping);
406 Finley_NodeMapping_free(in->Nodes->reducedDegreesOfFreedomMapping);
407 Paso_Distribution_free(in->Nodes->nodesDistribution);
408 Paso_Distribution_free(in->Nodes->reducedNodesDistribution);
409 Paso_Distribution_free(in->Nodes->degreesOfFreedomDistribution);
410 Paso_Distribution_free(in->Nodes->reducedDegreesOfFreedomDistribution);
411 Paso_Connector_free(in->Nodes->degreesOfFreedomConnector);
412 Paso_Connector_free(in->Nodes->reducedDegreesOfFreedomConnector);
413 in->Nodes->nodesMapping=NULL;
414 in->Nodes->reducedNodesMapping=NULL;
415 in->Nodes->degreesOfFreedomMapping=NULL;
416 in->Nodes->reducedDegreesOfFreedomMapping=NULL;
417 in->Nodes->nodesDistribution=NULL;
418 in->Nodes->reducedNodesDistribution=NULL;
419 in->Nodes->degreesOfFreedomDistribution=NULL;
420 in->Nodes->reducedDegreesOfFreedomDistribution=NULL;
421 in->Nodes->degreesOfFreedomConnector=NULL;
422 in->Nodes->reducedDegreesOfFreedomConnector=NULL;
423 }
424 }
425

  ViewVC Help
Powered by ViewVC 1.1.26